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Abstract 

Quantum field theory provides the framework for the most fundamental physical theories 
to be confirmed experimentally, and has enabled predictions of unprecedented precision. How- 
ever, calculations of physical observables often require great computational complexity and can 
generally be performed only when the interaction strength is weak. A full understanding of the 
foundations and rich consequences of quantum field theory remains an outstanding challenge. 
We develop a quantum algorithm to compute relativistic scattering amplitudes in massive </> 4 
theory in spacetime of four and fewer dimensions. The algorithm runs in a time that is poly- 
nomial in the number of particles, their energy, and the desired precision, and applies at both 
weak and strong coupling. Thus, it offers exponential speedup over existing classical methods 
at high precision or strong coupling. 
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1 Introduction 



Quantum field theory, by which quantum mechanics and special relativity are reconciled, applies 
quantum mechanics to functions of space and time — fields. In the quantum mechanics of point 
particles, dynamical variables become operators obeying suitable commutation relations, and the 
wavefunction encodes the probability that a particle is at position x. In quantum field theory, 
fields are promoted to operators satisfying commutation relations, and the wavefunction specifies 
the probability of each field configuration. 

Quantum field theory is the basis for the most fundamental physical theory to be established 
by experiment, namely, the Standard Model of particle physics. This quantum field theory en- 
compasses the strong, weak, and electromagnetic forces — all known forces except gravity — and 
explains all subatomic processes observed so far. Quantum field theory has also enabled predic- 
tions of unparalleled precision. For example, the most precise determination of the fine-structure 
constant (a measure of the strength of the electromagnetic force) uses quantum electrodynamics 
(QED) theory and the electron magnetic moment; it gives [1] 



in excellent agreement with the next most precise independent determinations [21 E]. 

To compare the predictions of quantum field theories with experimental results from particle 
accelerators, one typically calculates scattering amplitudes, using a perturbative expansion in pow- 
ers of the coupling constant, which quantifies the strength of interactions between particles. These 
perturbative calculations are greatly simplified through the use of Feynman diagrams. Neverthe- 
less, the number of Feynman diagrams needed to evaluate a scattering process involving n particles 
scales factorially with both the order in the perturbative expansion and n. Many-particle scattering 
processes can be phenomenologically important: for example, in a strongly coupled quantum field 
theory, the collision of two highly relativistic particles (with momentum much larger than mass, 
|p| ^> m) can produce a shower of n out ~ |p|/fra outgoing particles. 

In cases where the coupling constant is strong, such as at low energies in quantum chromody- 
namics (QCD), the quantum field theory governing nuclear physics, perturbative techniques break 
down. Strongly coupled quantum field theories can be studied on supercomputers via a method 
called lattice field theory. However, lattice field theory is useful only for computing static quantities 
such as mass ratios, and cannot predict dynamical quantities such as scattering amplitudes. 

In this work we develop a quantum algorithm to calculate scattering amplitudes in a massive 
quantum field theory with a quartic self-interaction, called (j) 4 theory. The complexity of our 
algorithm is polynomial in the time and volume being simulated, the number of external particles, 
and the energy scale. The algorithm enables the calculation of amplitudes to arbitrarily high 
precision, and at strong coupling. Specifically, the asymptotic scaling of the number of quantum 
gates, G W eak! required to sample from scattering probabilities in the weakly coupled, (d + 1)- 
dimensional theory with accuracy ±e is ( ^3.2p 





G we ak — * 



' /l\l-5+o(l) 
/ lx2.376+o(l) 



d = l 



d = 2 




1 \ 3.564+o(l) 



d = 3. 



For the strongly coupled theory, the asymptotic scaling is given in Table [T] of §3.21 
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There are at least three motivations for this work. First, there is a familiar motivation: the 
search for new problems that can be solved more efficiently by a quantum algorithm than by any 
existing classical algorithm. Secondly, one would like to know the physical limits of computation. 
What is the computational power of our universe? Since our most complete theory yet is the 
Standard Model, this is a question about the computational power of quantum field theories. 
Thirdly, we wish to learn more about the nature and foundations of quantum field theory itself. 
Indeed, history provides a precedent: Wilson discovered deep insights — of great conceptual and 
technical significance — about quantum field theory through his work on simulating it on classical 
computers [3J. 

The extensive efforts that are currently being made to build a quantum computer are largely 
motivated by the discovery of a small number of quantum algorithms that can solve certain com- 
putational problems in exponentially fewer steps (as a function of input size) than the best existing 
classical algorithms. The most famous example is Shor's quantum algorithm for prime factoriza- 
tion of an n-bit integer in 0(n 2 log n log log n) steps [5]. In comparison, the fastest known classical 
algorithm uses 2 ( nl/3 ( 1 ° sra ) 2/ ' 3 ) steps [6]. There has now been a large amount of effort devoted to 
finding additional quantum algorithms offering exponential speedup. This turns out to be very diffi- 
cult. At present, only relatively few problems have been discovered admitting exponential quantum 
speedup [7], including simulating nonrelativistic many-body quantum systems, estimating certain 
topological invariants such as Jones polynomials, and solving certain number-theory problems such 
as Pell's equation. 

The theory of classical computation was put on a firm foundation by proofs that Turing ma- 
chines, possibly supplemented by random-number generation, can efficiently simulate classical sys- 
tems, and, in particular, all other proposed models of universal classical computation. Thus, 
computational-complexity results obtained within the Turing machine model, or any of the equiva- 
lent models, are universally applicable, up to polynomial factors. The field of quantum computing 
originated with Feynman's observation that quantum-mechanical systems represent an apparent 
exception to this efficient mutual simulatability: simulating quantum-mechanical systems with n 
degrees of freedom on a classical computer requires time scaling exponentially in n [8]. Feynman 
conjectured that this exception was an indication of two classes of computation: classical and quan- 
tum, with all quantum models able to simulate each other efficiently, and all classical models able 
to simulate each other, but not quantum models, efficiently. 

In 1985, Deutsch defined quantum Turing machines and proposed them as a concrete model 
of a universal quantum computer able efficiently to simulate all others [9j- Subsequently proposed 
models of quantum computation, such as the quantum circuit model [lOj, continuous [11] and dis- 
crete [12] quantum walks, topological quantum computers [13j . and adiabatic quantum computers 
|14j . have all been proven to be polynomially equivalent to the quantum Turing machine model 
[I"5l [T6| [TTt [H2 [El [20]. Today, the model of quantum computation most widely used by theorists 
is the quantum circuit model. Feynman's original question about whether a universal quantum 
computer can efficiently simulate not only other models of quantum computation but also natu- 
rally occurring quantum systems has also been addressed: quantum circuits have been designed 
to simulate the evolution for time t of systems of n nonrelativistic particles using poly(n, t) gates 

[2U mi [sa M- 

Whether quantum circuits can efficiently simulate quantum field theories has remained an open 
question. Interactions between fields are nonlinear (that is, there are nonlinearities in the equations 
of motion), so the dynamics of an interacting theory is non-trivial to compute. Interactions are also 
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local, occurring at a point in spacetime, and hence there are infinitely many degrees of freedom 
even within a finite volume. Thus, known techniques for simulating the dynamics of many-body 
quantum- mechanical systems do not directly apply. If quantum field theories could not be simulated 
in polynomial time by quantum circuits, this would undermine the status of quantum circuits as 
a universal model of quantum computation and raise the possibility of quantum field computers, 
exponentially faster for certain tasks than even quantum circuits. 

Known classical methods for simulating quantum field theory require, on worst-case instances, 
exponential time to compute scattering amplitudes. Moreover, complexity theory provides reasons 
to believe that this is a fundamental barrier for all classical algorithms. It is generally believed 
that, in principle, a quantum computer could be built from known forms of matter. (The obstacles 
are mainly related to achieving high precision and good isolation from noise.) Thus, the quantum 
field theory governing known forms of matter (the Standard Model) should be sufficiently rich 
to model an idealized universal quantum computer. Hence, polynomial-time classical simulation 
of the Standard Model should be impossible, unless all quantum computations can be efficiently 
performed by classical computers. One could restate this claim in the language of complexity theory 
as a formal conjecture: the problem of simulating the Standard Model is BQP-hard, and therefore 
cannot be contained in P unless P = BQP. 

When one simulates classical physics, it is standard to discretize space using a lattice. The 
same thing can be done in quantum field theory, and this method underlies the lattice field theory 
calculations extensively used in supercomputer studies of quantum chromodynamics. To simulate a 
typical process at energy scale E, it is believed to suffice if one chooses a lattice spacing small com- 
pared with ^. However, discretizing a quantum field theory involves special difficulties unfamiliar 
from the classical context. A discretized system from classical physics converges to a continuum 
limit when one straightforwardly takes the lattice spacing to zero. In contrast, continuum limits 
of quantum field theories are achieved only through careful adjustment of the parameters of the 
Hamiltonian (or, equivalently, Lagrangian) as a function of the lattice spacing. This is the process 
of renormalization, which we discuss in £|2.1.2[ and it is an important consideration in the analysis 
of quantum (or classical) algorithms for simulating quantum field theories. In §5.31 we use methods 
of effective field theory to make a detailed analysis of the errors introduced to our simulation by 
discretization. 

Once discretized, a quantum field theory becomes essentially an ordinary many-body quantum- 
mechanical system, whose evolution can be efficiently simulated on quantum computers with estab- 
lished methods [25]. However, a full simulation includes, in addition to time evolution, preparation 
of a physically meaningful initial state, and final measurement of physically relevant observables. 
The state preparation and measurements depend strongly on the underlying physics, and must be 
analyzed on a case- by-case basis. In $3]we propose an adiabatic method for preparing wavepacket 
states of (f) A theory, and in §4.2l we analyze its complexity. We similarly construct observables mod- 
eling particle detectors by adiabatically continuing to the free theory ( §3.11 and !4.4p . Our analysis 
of state preparation and discretization are two of the most important technical contributions of 
this paper. 

No previous paper has addressed the quantum computation of scattering amplitudes or the 
convergence of quantum simulations to the continuum limit of a quantum field theory. The issue of 
gauge symmetries in qubit representations of lattice field theories has been discussed in [25]. There 
is also an extensive literature on how experimentally to construct Hamiltonians that approximate 
lattice gauge theories, in systems of atoms [551 [27] or superconducting qubits [28] . Atomic-physics- 
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based implementations have also been proposed for quantum field theories in curved spacetime 
[29j 1301 EU [32l [33] and superstrings [Ml [35]. Rather than addressing convergence to the contin- 
uum or preparation of wavepackets, these proposals focus on experimentally implementing lattice 
Hamiltonians, cooling to the ground state ("vacuum"), and measuring static quantities. 

The experimental implementations of Hamiltonians approximating lattice theories could be 
viewed as specialized quantum computers. However, they differ from the quantum circuit model 
considered here in that they are a form of analog, rather than digital, quantum computation. The 
proofs that digital quantum computers can in principle perform universal, scalable, fault-tolerant 
computation (e.g. [361 EH E] ) do not apply to analog quantum computers. The studies of analog 
and digital quantum simulation thus serve complementary purposes. By studying digital quan- 
tum simulation we may understand the ultimate abilities and limits of future large-scale quantum 
computers, whereas analog quantum simulation offers more specialized and less scalable techniques 
that are implementable with present-day technology. 

Some of the techniques discussed in this paper may be useful for designing quantum algorithms 
to do things other than the simulation of high-energy physics. In particular, suppose one wishes to 
simulate the quantum dynamics of a crystal lattice. One could do this directly using the standard 
quantum algorithms for quantum simulation, and the resources would scale polynomially in the 
number of lattice sites. However, if the phenomena being studied are characterized by much longer 
length scales than the crystal lattice spacing, this direct approach may be unnecessarily costly. 
Instead, one could use the renormalization group to obtain an effective Hamiltonian on a coarse- 
grained lattice, and then simulate that using the standard techniques. 

2 Background 

Our presentation is intended to be self-contained and accessible to a broad audienc^ll Section 12.11 
describes some basic concepts of quantum field theory drawn upon in this work. For simplicity, the 
discussion is restricted to the scalar (ft 4 theory, which is the focus of this paper. Section 12.21 gives 
an overview of the standard techniques in quantum computing that are used in our algorithm. For 
reference, notation introduced in this and subsequent sections is tabulated in Appendix A. 

2.1 Quantum Field Theory 

The scalar (ft 4 quantum field theory in the continuum limit is relevant to diverse physical phenomena. 
In the Standard Model of particle physics, the Higgs boson endows mass to other particles. Its self- 
interactions are described by the (four-component) (ft theory. Additionally, the Euclidean versions 
of (ft 4 theories have been highly successful in describing critical phenomena and universality in 
statistical-mechanical systems. 

The simplicity of (ft theory makes it suitable for studying challenging formal aspects of quantum 
field theory. The existence of the (ft 4 continuum theory has been rigorously established in two 
[39l l40l I4T1 [42j [43] and three [HI 05] (spacetime) dimensions. The perturbative expansion for 
the scattering matrix has been shown to be asymptotic in two |46[ |4"T] and three |48j dimensions. 
Furthermore, in those dimensions the Schwinger functions (Euclidean correlation functions) can be 
reconstructed from their perturbative expansions, that is, they are Borel summable [491 150[ 1ST]. In 
contrast, in five or more dimensions, (ft 4 theory is trivial [521 153j. that is, the continuum limit is 

For a concise description of the main ideas, see |38j . 



6 



equivalent to a free theory. In the borderline case of four dimensions, obtaining a rigorous proof is 
still an open question, but the theory is believed to be trivial. Despite this, it is still an important 
theory: it is simply regarded as an effective field theory valid below some high-energy cutoff. 

One possibility for attempting to construct a continuum theory in an infinite volume is to define 
a theory on a finite lattice and then study the limits as the lattice spacing is taken to zero and the 
volume to infinity. In the following, we present a lattice theory and its continuum limit in a manner 
directly corresponding to the basis of our algorithm. From a field theorist's perspective, the setup 
is therefore unorthodox in several respects. Only the spatial dimensions of the (Minkowski-space) 
theory are put on a lattice, as opposed to all spacetime dimensions of a Euclidean theory. The 
Schrodinger picture is used rather than the interaction picture, and the Hamiltonian formalism 
rather than the Lagrangian formalism. As is conventional in quantum field theory, we use units 
where h = c = 1, so that all quantities have units of some power of mass. At the end, any necessary 
factors of H and c can be reinserted by dimensional analysis. 

2.1.1 Lattice </> 4 Theory 

Let be a cubic lattice in d spatial dimensions with length L, volume V = L d , lattice spacing 
a, and V = {L/a) d lattice sites. Also let L = L/a. For simplicity, we assume periodic boundary 
conditions, that is, 

n = ol\. (2) 

For each x £ 0, let </>(x) be a continuous, real degree of freedom interpreted as the field at x, 
and 7r(x) the canonically conjugate variable to </>(x). By canonical quantization, these degrees of 
freedom are promoted to Hermitian operators with the following commutation relations: 

[0(x),vr(y)] = ia- d 5 x , y l, 
[0(x),0(y)] = [vr(x),vr(y)] = 0. (3) 

Here, the factor of a~ d provides dimensions of (length) - ^, matching a Dirac delta function in the 
(i-dimensional continuum. (J) 4 theory on the lattice f2 is defined by the following Hamiltonian: 



H=J2a d 



ivr(x) 2 + \ (V a cf>) 2 (x) + \mlm 2 + ^0(x) 4 



(4) 



V a (f> denotes a discretized derivative, so that 

(v.«'w-£ ( fetisl^w V, (5) 

j=i 

where r\, . . . , are unit vectors such that afi, . . . , af^ form the basis of the lattice The ground 
state of H is called the vacuum and denoted |vac). It can be interpreted as the state in which there 
are no particles. The covariance matrix 

G(x- y ) = (vac|0(x)0( y )|vac) (6) 

is very useful in quantum field theory, and is referred to as the equal-time propagator. In the free 
theory (Ao = 0), the equal-time propagator can be exactly evaluated. Let T be the momentum-space 
lattice corresponding to fl: 

r = ^z d L . (7) 
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Then, 



G^x-y) 



1 1 



^ L d 2w(p) 
per yt ^ J 



Jp-(Xi-Xj) 



(8) 



Note that J2xeQ a<i ~ * I ^ x asffi_> "i an d S p er Ja ~ * I pLy • as L — )• oo. 

Let be the Hamiltonian of with Aq = 0. One can exactly solve for the spectrum of 



by rewriting it in terms of creation and annihilation operators, a P and a p . For each pGf, let 



a d e~ ip - 



" (P) -X) + Zj^vr(x) 



where 



One can verify that 



w(p) 



l + n 2^2'' 



2 ap-r 

j=l 



x 



7r(x) 
H (0) 

Op, Oq 

[a p ,a q ] 




per 

^^(p)a p o p + ^°)l 



per 

-d 



0, 



0) 
(10) 

(11) 
(12) 
(13) 

(14) 
(15) 



with = X] p6 r ^ w (p)- In the limit 

a — > we have [a p ,a q ] — > (27r) rf (5^^ (p — q), which is the 
standard relation for creation and annihilation operators in the continuum. 

Let |vac(0)) denote the ground state of F (0) ("the free vacuum"). By ([13]), CHD, and (fT5|) . 
one sees that L _d//2 a P |vac(0)) is an eigenstate of of unit normal with energy w(p). This state 
can be interpreted as containing a single particle of momentum p. Applying additional creation 
operators with sharply defined momentum yields an eigenstate of H^ ' whose energy is simply 
the sum of the energies of the particles created. In this way, one obtains the entire spectrum of 
this non-interacting theory. A perfectly momentum-resolving particle detector can be modeled 
by the observable L~ d a P a p , whose eigenvalues are integers that count how many particles are in 
momentum mode p. 

In accordance with the uncertainty principle, L _d/ ' 2 a P |vac(0)) represents a particle completely 
delocalized in space. The state 0(x)|vac(O)) is interpreted as a single particle localized at x (up to 
normalization). Because a p |vac(0)) = 0, 0(x)|vac(O)) = a£|vac(0)), where 



1 



-?px 



1 



^ L d 
per 



2w(p) 



(16) 



The standard relativistic normalization convention is not convenient in our context. 
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In the nonrelativistic limit rriQ 3> |p|, the factor of in f|16|) becomes approximately constant, and 
one recovers the familiar nonrelativistic Fourier relation between position-space and momentum- 
space wavef unctions. 

We are able to solve for the spectrum of in terms of non-interacting particles with sharply 
defined momentum. In an interacting quantum field theory, the momentum of individual particles 
is no longer conserved. Nevertheless, total momentum is conserved, and thus the single-particle 
subspace can be analyzed similarly to the non-interacting case. Specifically, by starting with the 
single-particle momentum-p eigenstate of the non-interacting theory, L -<i / 2 ap|vac(0)), and then 
adiabatically turning on Ao, one obtains an eigenstate of H, which can be interpreted as a single 
particle of momentum p of the interacting theory. 

2.1.2 Continuum Limit 

In the preceding section we described a lattice field theory. The quantum mechanics of such a 
system is mathematically well defined. However, a lattice theory lacks translational, rotational, 
and Lorentz invariance. Thus one is led to try to construct a continuum limit of a sequence of 
lattice theories with successively finer spacing, much as one obtains integrals by a sequence of more 
finely discretized Riemann sums. 

The naive attempt at constructing a continuum limit in which one simply takes a — > in the 
definition of the lattice Hamiltonian fails to yield convergent answers to physical questions. Instead, 
one should consider a sequence of Hamiltonians of the form (j4]) on successively finer lattices, where 
uiq and Ao are functions of a. One aims to show that, for a suitable choice of the a dependence 
of rriQ and Ao, the entire theory converges to some meaningful limit as o — > 0. This is known as 
renormalization. 

Although the existence of a continuum limit of a sequence of (Euclidean) <\>\ 3 lattice theories 
has been shown rigorously in [54j . all that has been demonstrated for most physically interesting 
quantum field theories is perturbative renormalizability, namely, that physical quantities are finite 
to all orders in perturbation theory. Physical (renormalized) quantities are then calculated as a 
perturbative series in the coupling. In our analysis of weak coupling, we use expressions obtained 
from perturbation theory. Specifically, in §5.11 we calculate the physical mass m up to second 
order in Ao- Inverting this expression and the analogous one for A yields prescriptions for choosing 
the bare parameters mo and Ao as functions of a to achieve given physical parameters A and m 
(corresponding to what is measured by experiment). 

For the purpose of simulation, in particular the analysis of algorithmic complexity, one must 
ask not only whether the sequence of Hamiltonians converges to a continuum limit, but also how 
quickly it converges to that limit. Such questions are addressed in §5.31 Similarly, the rate of 
convergence to an infinite-volume limit is studied in §5.41 

2.1.3 Phase Transitions 

Consider the continuum c/> 4 theory in the infinite- volume limit in D = 2 or 3 spacetime dimensions. 
At certain values of the parameters ttt-q and Ao, the particle mass m vanishes. These points form a 
critical curve, across which a quantum (zero-temperature) phase transition occurs, corresponding 
to the spontaneous breaking of the (j) — > — (f> symmetry of the theory. 

The existence of such a phase transition was shown rigorously in [551 l56l [57] . As the system 
approaches it, thermodynamic functions and correlation functions exhibit power-law behavior, as 
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is characteristic of a second-order phase transition. In particular, for constant rriQ 



m ~ | A - A c 



(17) 



where A c , the critical value of the coupling, depends on jtIq. The critical exponent v depends only 
on the dimension and is bounded below: v > 1/2 |55|, I57j. 

Empirically, it has been found that systems with second-order phase transitions can be classified 
into universality classes. Within each class, critical exponents are universal, taking the same values 
for all systems, even though the systems may be vastly different in their microscopic interactions. 
(That the critical behavior depends only on the symmetries and the spatial dimensionality of the 
Hamiltonian is explained by the concept of the renormalization group.) The (ft 4 theory is believed 
to be in the same universality class as the Ising model, for which 



The value above for D = 3 has also been obtained directly in the (ft 4 theory by Borel resummation 

Pi- 
rn general, the non-perturbative regime, in which the coupling is sufficiently strong that pertur- 
bation theory fails (and hence also called the strong-coupling regime), is expected to occur in the 
vicinity of the phase transition, with the inter-particle force maximum at the critical value of the 
coupling. Correspondingly, the dimensionless ratio \/m 4 ~ D (D = 2, 3) will become large: indeed, 
this is consistent with the mass shrinking to zero and the coupling approaching a nontrivial infrared 
fixed point. 

In D = 4 dimensions, in contrast, the believed triviality of the continuum (ft 4 theory implies that 
there is no non-trivial fixed point of the renormalization group and hence no phase transition as 
one varies (toq, Ao). Moreover, triviality places bounds on the maximum value of the renormalized 
coupling [59|. In particular, strong coupling requires pa to be 0(1): in the continuum-like regime, 
renormalized perturbation theory should be valid. 

2.2 Quantum Computing 

The quantum circuit model is a convenient framework for describing quantum computation]! Quan- 
tum circuits are in many respects analogous to classical circuits. In classical circuits, logic gates 
operate on bits, each in the state or 1. There exist universal sets of gates (on their own able 
to implement any Boolean function): examples are {AND, OR, NOT} and {NAND}. Similarly, in 
quantum circuits, quantum gates operate unitarily on qubits, each of which is in a linear superpo- 
sition of two basis states (|0) and |1)). Arbitrary single-qubit and controlled-NOT gates together 
can perform any unitary operation. Furthermore, there exist smal@ sets of gates that are universal, 
in the sense that they can approximate any unitary operation to arbitrary accuracy. 

In the following subsections we review the standard primitives from the theory of quantum 
circuits that are used in our algorithm. The number of gates needed to implement a unitary 
transformation is used as a measure of its running time, although the actual running time might 
be much shorter if the computation can be highly parallelized. 

3 An introduction to the basic notions of quantum circuits can be found in [60] . 

There exist two-qubit gates that are universal by themselves, but the standard universal gate set consists of three 
gates: controlled-NOT, Hadamard, and n/8 60 . 




D = 2, 
D = 3. 



(18) 
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2.2.1 Quantum Fourier Transform 



For any / : {0, 1, 2, . . . , JV - 1} ->■ C, the function / : {0, 1, . . . , N - 1} ->■ C given by 

N-l 

/» = -7W E e_27riifc/iV /(i) (19) 

is called the discrete Fourier transform of /. The linear transformation / — > f on the ^-dimensional 
complex vector space of all functions / : {0, 1, 2, . . . , N — 1} — > C is unitary. It can be implemented 
on the amplitudes of an arbitrary state with O (log 2 N log log N log log log N) quantum gates [5]. 



2.2.2 Reversible Circuits and Phase Kickback 

Because the quantum time evolution of a closed system is unitary, information cannot be erased in 
a quantum circuit. Suppose we are given a classical circuit computing some function / : {0, l} n — > 
{0, 1, 2, . . . , M — 1}. (We think of the output as an integer, although it is of course represented as a 
string of m = |~log 2 M] bits, called a "register".) Unless / is injective, there is no unitary operator 
Vf such that Vf\x) = \f(x)}. However, one can always define a unitary operator Uf such that, for 
any y£ {0,1,2,..., M - 1}, 

U f \x)\y) = \x)\{y + f{x)) mod M) . (20) 

Furthermore, the number of quantum gates needed to implement this unitary operator can never 
exceed the number of classical gates needed to compute / by more than a constant factor [61] ED] ■ 
If we initialize the output qubits to the zero state, then applying Uf yields f(x) written into the 
values of the bits: 

U f \x)\0) = \x)\f(x)). (21) 
If we instead prepare the output qubits in the state 

M-l 

\RM) = ^=Y, e ~ l27TV/M \y^ ( 22 ) 

v M 



then applying Uf "kicks back" f(x) into the phase that is, 

U f \x)\R M ) = ^^l M \x)\R M ). (23) 
The state \Rm) can be prepared efficiently by means of a quantum Fourier transform. 



2.2.3 Phase Estimation 

Consider a quantum circuit implementing a unitary transformation U on n qubits, and an eigenstate 
\9) such that U\6) = e lS \6). One can use quantum Fourier transforms to measure 6 to m bits of 
precision using 0(m 2 log m log log m + T) quantum gates, where T is the number of gates needed 
to implement U 2 ™ [63]. Applying the phase estimation circuit to an arbitrary state yields an 
approximate measurement in the eigenbasis of U. Eigenvalues closer together than ~ 2 _m are not 
distinguished. In many physical applications, one wishes to measure in the eigenbasis of a Hermitian 
operator Q. The eigenbasis of Q is the same as the eigenbasis of the unitary operator e 1 ® 1 . Thus, 
in cases where one can implement the unitary transformation e 1 ® 1 by an efficient quantum circuit, 
the problem of measuring Q reduces to the standard phase-estimation problem. More discussion 
of the implementation of e 1 ® 1 is given in the next section. 
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2.2.4 Simulating Hamiltonian Time Evolution 



Given a Hamiltonian, we wish to find a quantum circuit of few gates implementing the corresponding 
unitary time evolution. The details of how to do this depend on the specifics of the Hamiltonian. 
Thus for concreteness we specifically consider the Hamiltonian H defined in Q. The simulation 
method described here was introduced in [211 [23] . 

To simulate time evolution according to H, we decompose H as 

H = H^ + H^ (24) 

H * = jE A (xf, (25) 



2 



^ = E 

xesi 



a d 



\ (V a 0) 2 (x) + ^(x) 2 + ^>(x) 4 



(26) 



acts on the computational basis states simply by inducing a phase. This phase (namely, 
Sxen ad \ a4>) 2 ( x ) + ~2 L 4 I ( X } 2 + ^o</ ) ( x ) 4 ) is not hard to compute; on a classical computer one 
could compute it using 0(V) gates. Thus, using the method of phase kickback, one can implement 



H^St 



for any 5t using 0(V) quantum gates (j j2.2.2j) . 



Similarly, we can simulate e~ lH ^ St for any 5t by first Fourier transforming each of the V registers. 
(Because each register contains only logarithmically many qubits ( §4.ip . these Fourier transforms 
use logarithmically many quantum gates each, and thus they make only a small contribution to 
the overall complexity of the algorithm.) This brings the state into the eigenbasis of the com- 
plete set of commuting observables {7r(x)|x 6 (]}. In this basis, e~ lHnSt simply induces a phase 
Sxen o rf 7r(x) 2 (5t. Computing this phase has complexity O(V), and therefore phase kickback imple- 
ments e~ lH ^ St using 0(V) quantum gates. Afterwards, inverse Fourier transforms bring the system 
back into the computational basis, which is the eigenbasis of {</>(x)|x G O}. 

Given the ability to simulate e - *^ 5 * and e~ lH7vSt , one can simulate e - l { H <t> +H ^) t using Trotter's 
formula [24] . 

e -i(H++H«)t = ( e -iH^ln e -iH^tln\ n + Q(t 2 [H n , H^}/n) . (27) 



Thus, to achieve error esT one alternatingly performs time evolutions of and H n with n 
t 2 \\ [H w , Hf/,] ||/esT steps. Each step requires 0(V) gates. = 0(V) because of locality. 

Thus, the total complexity is 0(t 2 V 2 ). 

Fortunately, one can systematically construct higher-order Suzuki- Trotter formulae that yield 
better scaling in t. Using a /c th -order Suzuki- Trotter formula, one can simulate time evolution 
for time t on a lattice of V sites using a quantum circuit of 0((tV) 1+ 2fc) gates ( §4.3)) . Note that 
our analysis of the V scaling of quantum circuits to simulate spatially local Hamiltonians on large 
lattices is, to our knowledge, novel. This may be of independent interest for the simulation of lattice 
Hamiltonians arising, for example, in condensed-matter physics. 



2.2.5 Adiabatic State Preparation 

Preparing an arbitrary state on n qubits requires exponentially many quantum gates Even 
the restricted problem of preparing ground states of local Hamiltonians is, for worst-case instances, 
intractable for quantum computers. (More precisely, this problem is QMA-complete |65) 166 ) [67]. ) 
However, for certain Hamiltonians one can solve the problem of preparing the ground state using 
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polynomially many gates by simulating adiabatic time evolution, that is, time evolution according 
to a Hamiltonian with slowly time- varying parameters. 

Let .ffinit be a Hamiltonian whose ground state is easy to prepare. Suppose there exists some 
Hamiltonian H(s) such that H(0) = flhut and H(l) = -£ffi na i, and let \ip(s)) be the ground state 
of H(s). According to the adiabatic theorem, by starting with \ip(0)), and evolving according 
to the time-dependent Hamiltonian H(t/r) for sufficiently large r, one obtains at time r a state 
approximately equal to Let 7 be the energy gap between the ground state and first excited 

state. Quantitative versions of the adiabatic theorem have been proven that show it always suffices 
to choose t = 0(l/7 3 ) [68l [69] , and in certain cases it suffices to choose r = 0(l/7 2 ) [70J. 

For efficient state preparation within the quantum circuit model, it suffices if 7(5) is not too 
small and time evolution governed by H(s) can be simulated efficiently. 

3 Quantum Algorithm 

The set of field operators {</>(x)|x G Q} forms a complete set of commuting observables. We 
represent the state of our lattice field theory by devoting one register of qubits to store the value of 
the field at each lattice point. One can simulate a scattering process at energy scale E with fidelity 
1 — e to the exact state, with each register consisting of O (log (^jr)) qubits ( §4.ip . 
Schematically, our quantum algorithm proceeds as follows. 

1. Use the method of Kitaev and Webb [71j to prepare the ground state of the free theory 
(A = 0). 

2. Excite wavepackets of the free theory. 

3. Evolve for a time r during which the interaction is adiabatically turned on. This yields 
wavepackets of the interacting theory. 

4. Evolve for a time t in which scattering occurs. 

5. Evolve for a time r during which the interaction is adiabatically turned off. This brings us 
back to the free theory, where interpreting field states in terms of particles is straightforward. 

6. Measure the number operators of the momentum modes of the free theory. 
3.1 Description and Discussion 

We next discuss each step of the quantum algorithm in more detail. The full analysis of their 
complexity scaling occupies $H 

Instead of {0(x)|x G Q}, one could use {apO p |p G T} as the complete set of commuting 
observables. In this case, the qubits store the occupation numbers of momentum modes. However, 
simulations in the field representation seem to be more efficient than simulations in the occupation- 
number representation. The primary reason for this is that, upon expanding the c/> 4 operator 
in terms of creation and annihilation operators, one obtains an expression that is nonlocal in 
momentum space. This makes the simulation of eT %Ht using Suzuki- Trotter formulae somewhat 
inefficient (although still polynomial-time). Therefore, throughout this paper we consider only 
simulations based on the field representation. 
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1. Improving upon the efficiency of earlier, more general, state-construction methods |2H I72j. 
Kitaev and Webb developed a quantum algorithm for constructing multivariate Gaussian 
superpositions [71]. The main idea of this algorithm is to prepare a V-dimensional multivariate 
Gaussian wavefunction with a diagonal covariance matrix, and then reversibly change basis to 
obtain the desired covariance matrix. For large V, the dominant cost in Kitaev and Webb's 
method is the computation of the LDL T decomposition of the inverse covariance matrix, 
where L is a unit lower triangular matrix, and D is a diagonal matrix. This can be done 
in 0(V 2 ' 376 ) time with established classical methods [731 [74]. (The notation f(n) = 0(g{n)) 
means f(n) = 0(g(n) log c (n)) for some constant c.) The computation of the matrix elements 
of the covariance matrix itself is easy because, for large V, the sum (jSJ) is well approximated 
by an easily evaluated integral. 

2. The expression 0(x)|vac) for a single-particle state is not directly useful for quantum com- 
puting because </>(x) is not unitary. Instead, given a position-space wavefunction ij) (with 
normalization X^xeo ^I^MI 2 = 1 so that it matches the dimensions of a continuum wave- 
function) let 

h/>>=4,|vac(0)>, (28) 

where 

a{=#)^«^(x)4, (29) 

4 is defined in (|16|) . and T}{$) is the normalization factor that ensures o^a^|vac(0)) = |vac(0)). 
Introduce one ancillary qubit, and let 

fl^ = a^®|l)(0| + tty®|0)(l|. (30) 

One can easily verify that the span of |vac(0))|0) and is an invariant subspace, on which 

H^jj acts as 

ify|vac(0)>|0> = (31) 
H^)\l) = |vac(0))|0>. (32) 

Thus 

e-^ 7r/2 |vac(0))|0) = (33) 

Hence, by simulating a time evolution according to the Hamiltonian H^, we obtain the desired 
wavepacket state \tp), up to an irrelevant global phase and extra qubit, which can be discarded. 
After rewriting Hj, in terms of the operators 0(x) and vr(x), one sees that simulating is 
a very similar task to simulating H, and can be done with the same techniques described in 
and EM 

For a spatially localized wavepacket, that is, ifi with bounded support, is a quasilocal 
operator ( §4.5f) : it can be written in the form ^ x [/(x)7r(x) + g(x)0(x)], where / and g 
decay exponentially with characteristic decay length 1/mo outside the support of ip. One can 
thus choose some c\ 3> 1 and set / and g to zero outside a distance ci/mo from the support 
of ip. The resulting operator will be fully local and an exponentially good approximation 
to H^. The time evolution according to this local approximation to can be simulated 
with complexity independent of V. Furthermore, this quasilocality shows that, to create 
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wavepackets of additional particles, we can simply repeat this procedure with different, well 
separated, choices of the position-space wavefunction t/j. The only errors introduced at this 
step are due to the finite separation distance 5 between wavepackets, and are of order ei oc ~ 
e ~5/m^ (However, our wavepackets have a constant spread in momentum, and thus differ from 
the idealization of particles with precisely defined momenta.) The wavepacket preparation 
thus has complexity scaling linearly with n; n , the number of particles being prepared, and 
necessitates a dependence V ~ n- m log(l/ei oc ). 

3. To obtain a wavepacket of the interacting theory, we start with the wavepacket of the free 
theory, constructed in the previous step, and then simulate adiabatic turn-on of the interac- 
tion. An adiabatic process of time r can be simulated by a quantum circuit of 0((rV) 1_ 
gates implementing a fe th -order Suzuki- Trotter formula (£ j4.3p . 

For < s < 1, let 



' 2k 



ivr(x) 2 + ^(V a 0) 2 (x) + ±mg( a )0(x) 2 + ^ 



x 



(34) 



with Ao(0) = 0. If we started with an eigenstate of the free theory H (0) and applied the time- 
dependent Hamiltonian H(t/r) for time r then, by the adiabatic theorem, we would obtain a 
good approximation to the corresponding eigenstate of H(l), provided r was sufficiently large. 
However, a wavepacket is a superposition of eigenstates with different energies. These acquire 
different phases during the adiabatic state preparation. Physically, this means the wavepack- 
ets propagate and broaden. Propagation during adiabatic state preparation is undesirable 
because the particles being prepared could collide and scatter prematurely. Broadening of 
wavepackets decreases the efficiency of the algorithm because very diffuse wavepackets largely 
pass through each other without significant scattering. (For diffuse wavepackets, the expected 
distance between particles is large even when the wavepackets are on top of each other.) In this 
case, many repetitions of the simulation are necessary before observing interesting scattering 
events. 

We can correct the dynamical phases by interspersing the adiabatic state preparation with 
backward time evolutions, thereby counteracting the propagation and broadening of wavepack- 
ets. Specifically, we follow an adiabatic path H (s) from s = to s = 1. To undo the dynamical 
phase, we divide the total adiabatic evolution into J steps. Before and after each step we 
apply time-independent Hamiltonians as follows. For j = 0, 1, . . . , J — 1, let 



Mi 



exp 



iH 



J 23 



Uj exp 



J 2 J 







Uj = T | exp 


—i / H{s)rds 







(35) 
(36) 



where T{-} indicates the time-ordered product. (Uj is the unitary time evolution induced by 
H(t/r) from t = to t = ■) The full state-preparation process is given by the unitary 



J-i 

M=J[M j . 

3=0 



(37) 
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We suppress the dynamical phases by choosing J to be sufficiently large. The choice of a 
suitable "path" Ao(s), mg(s), and the complexity of this state-preparation process depends in 
a complicated manner on the parameters in H (£ )4,2p , 

4. The time evolution e - iH {i)t can be implemented with O (^(tV) 1+ ^^j gates via a /c th -order 
Suzuki- Trotter formula ( §4.3p . 

5. The adiabatic turn-off of the coupling is simply the time-reversed version of the adiabatic 
turn-on. The complexity is no higher than that of the adiabatic state preparation^!. 

6. By the method of phase estimation, measurement of the number operator L^ d a\ ) a v> reduces 
to the problem of simulating e tL da p a P t for various t. Thus, for a given p, L~ d a\,a p can 
be measured with O (v 2+ 2k^j quantum gates via a A; th -order Suzuki- Trotter formula ( §4.4j) . 
Furthermore, if we instead simulate localized detectors, the computational cost becomes in- 
dependent of V (much as the computational cost of creating local wavepackets is independent 
of V ( §4.5|) ). but the momentum resolution becomes lower, as dictated by the uncertainty 
principle. 

The allowable rate of adiabatic increase of the coupling constant during state preparation is 
determined by the physical mass of the theory. In the weakly coupled case, this can be calculated 
perturbatively, as is done in £ 15,1.21 In the strongly coupled case, it is not known how to do the 
calculation. Thus one is left with the problem of determining how fast one can perform the adiabatic 
state preparation without introducing errors. Fortunately, one can easily calculate the mass using 
a quantum computer, as follows. First, one adiabatically prepares the interacting vacuum state 
at some small Ao, and measures the energy of the vacuum using phase estimation. The speed at 
which to increase Ao can be chosen perturbatively for this small value of Ao- Next, one adiabatically 
prepares the state with a single zero-momentum particle at the same value of Ao, and measures 
its energy using phase estimation. Taking the difference of these values yields the physical mass. 
This value of the physical mass provides guidance as to the speed of adiabatic increase of the 
coupling to reach a slightly higher Ao- Repeating this process for successively higher Ao allows one 
to reach strong coupling, while always having an estimate of mass by which to choose a safe speed 
for adiabatic state preparation. In addition, mapping out the physical mass as a function of bare 
parameters (hence, for example, mapping out the phase diagram) may be of independent interest. 

3.2 Efficiency 

In complexity theory, the efficiency of an algorithm is judged by how its computational demands 
scale with the problem size or some other quantity associated with the problem's intrinsic difficulty. 
An algorithm with polynomial-time asymptotic scaling is considered to be feasible, whereas one 
with super-polynomial (typically, exponential) scaling is considered infeasible. This classification 
has proved to be a very useful guide in practice. The results stated below can be roughly summarized 
as follows: the calculation of quantum field-theoretical scattering amplitudes at high precision or 
strong coupling is infeasible on classical computers but feasible on quantum computers. 

5 Readers of i]4.2l may notice that, in the strongly coupled case, the matrix element for particle splitting, which 
appears in the numerator of the diabatic error, gets multiplied by n out . However, this is more than compensated by 
the smaller momentum and hence larger energy gap against splitting. 
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Traditional calculations of QFT scattering amplitudes rely upon perturbation theory, namely, 
a series expansion in powers of the coupling (the coefficient of the interaction term), which is taken 
to be small. A powerful and intuitive way of organizing this perturbative expansion is through 
Feynman diagrams, in which the number of loops is associated with the power of the coupling. A 
reasonable measure of the computational complexity of perturbative calculations is therefore the 
number of Feynman diagrams involved. The number of diagrams is determined by combinatorics, 
and grows factorially with the number of loops and the number of external particles. If the coupling 
constant is insufficiently small, the perturbation series does not yield correct results. There are then 
no feasible classical methods for calculating scattering amplitudes, although lattice field theory can 
be used to obtain static quantities, such as mass ratios. For computing the time evolution of the 
state vector on a classical computer nonperturbatively, no better method is known than discretizing 
the field and solving the Schrodinger equation numerically. This method is infeasible because the 
dimension of the Hilbert space is exponentially large. 

Even at weak coupling, the perturbation series is not convergent, but only asymptotic: the error 
in the A rth -order sum Ylk=l ^k9 k satisfies 

N 

N+l 



M-J2b k g 



k=l 

Since the coefficients grow as (75J EH] 



0( 5 JV+i ) as g^O. (38) 



b k ~ k\c k k C2 , (39) 



for some constants c\ and C2, there is a maximum possible precision, corresponding to truncation 
of the series around the (l/g , ) th term. 

Consequently, our quantum algorithm should have an advantage in the non-perturbative regime, 
or if very high precision is required. Thus, we analyze the asymptotic scaling of our quantum algo- 
rithm when simulating weakly coupled theories at arbitrarily high precision and strongly coupled 
theories arbitrarily close to the quantum phase transition. In the strongly coupled case, we also 
consider the scaling as a function of the momenta of the incoming particles. As the energy of 
the incoming particles becomes larger, the maximum number of kinematically allowed outgoing 
particles correspondingly increases, thereby making the problem potentially more computationally 
difficult. In the weakly coupled case, processes producing large numbers of outgoing particles are 
suppressed even at high energy, because they arise only at high order in perturbation theory. 

In the weakly coupled case, we wish to determine the complexity of our algorithm as a function 
of precision. We quantify this by demanding that any standard physical quantity a extracted from 
the simulation (for example, a scattering cross section) satisfy 

(1 - e)(Texact < cr < (1 + e)<7 exact . (40) 

To analyze the scaling of our algorithm with e, we first consider errors due to spatial dis- 
cretization. The effect of spatial discretization is captured by (infinitely many) additional terms 
in the effective Hamiltonian ( §5.3p . Truncation of these terms alters the calculated probability of 
scattering events. In particular, the two dominant extra terms in the effective Hamiltonian are 
Yliddf'fi — 0^x^ an d <^ 6 terms, arising from discretization of (V a </>) 2 and quantum effects, respec- 
tively. The coefficient of the 4>d^(j) term is 0(a 2 ), and the coefficient of the </> 6 term is 0(a 5 ~ d ), 
so that the former dominates for d = 1,2, whereas the latter makes a comparable contribution for 
d = 3. (To improve the scaling, one can use better finite differences to approximate the derivative, 
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and/or include the </> 6 operator. However, renormalization and mixing of the coefficients make this 
idea more complicated than it is in standard numerical analysis.) 

We now describe the errors induced by neglecting the effective (j)d^(p term. The analysis of </> 6 
errors is analogous. Let a be a scattering cross section induced by the Hamiltonian H of ([4]), and 
let a' be the corresponding scattering cross section induced by 

H> = H -£a rf !^>(x) . (41) 

By standard arguments from effective field theory ( §5.2fJ5T3j) . 

a' = a (l + cf + 0(c 2 )), (42) 

where / is a function of the momenta and masses of the particles involved in the scattering process, 
c = 0{a 2 ) and / is independent of a, so the errors induced by neglecting (j)d^(j) ar e of order a 2 . 
Similarly, the errors induced by neglecting the effective <^> 6 term are of order a 2 or smaller for 
D = 2,3, 4. 

For the total error to satisfy (|40p . each individual source of error must be at most O(e). Thus, to 
ensure the spatial discretization errors are sufficiently small, we set a ~ y/e. This choice of a affects 
the complexity of the preparation of the free vacuum, and the complexity of the Suzuki- Trotter 
time evolutions. Because V = the Kitaev-Webb state preparation uses 

G prep = 0(V 2 ' 376 ) = 0(a- 2 ' 37M ) = 0(e- L18M ) (43) 

quantum gates. Among the various Suzuki- Trotter time evolutions in our algorithm, the most 
time-consuming is the adiabatic transition from the free theory to the interacting theory. We thus 
substitute a ~ \fk into (|124p and find that implementing this process with a A; th -order Suzuki- 
Trotter formula usea^l 

G adiaba tic ~ [f^^J) ^ > ^ = 1, 2, 3 (44) 

quantum gates, where, by definition, the adiabatically produced state has an inner product of at 
least 1 — e a( j with the exact state ( §4.2p . It is efficient to use Suzuki- Trotter formulae with large k, 
in which case the terms in the exponents such as Jfc become very small. To simplify our exposition, 
we henceforth use the standard "little-o" notation, defined as follows: 

/(n) = o{g{n)) if and only if lim f(n)/g(n) = 0. (45) 

(In the case of e scaling it is of course 1/e that is taken to infinity in the asymptotic notations.) In 
this language, (j4"4"|) becomes 

Gadiabatic ~ (e^Vd) ^ ^ . d = 1, 2, 3. (46) 

To obtain the full scaling with e of the adiabatic state preparation, we must next determine 
the relationship between e and e a d- For small e a d) the adiabatically prepared and exact states 
yield nearly identical probability distributions for all possible measurements. (More precisely, the 



Whether we use (|124|l or (|125[) affects only the scaling with V. 
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total variation distance from the exact probability distribution is at most 0(e a d).) Thus, setting 
£ad = 0{e) is certainly sufficient to satisfy (j4"0j) . We then have 



G. 



adiabatic 



-l-d/2+o(l) 



(47) 



Comparing ([4*3]) with ([IT)) , one sees that in d = 1 the adiabatic state preparation is the dom- 
inant cost, whereas in d = 2,3 the preparation of the free vacuum dominates. This leaves a final 
asymptotic scaling of 



G 'total — 0(G a( jiabatic + Gp Iep ) — < 



(1) L5+0(1) , d=l, 

(i) 2 - 376+o(1) , d = 2, 



(48) 



So far, in determining the e scaling of our algorithm, we have considered only errors due to 
spatial discretization and errors due to imperfect adiabaticity. We now argue that the remaining 
sources of error make negligible contributions to the overall e scaling, which are already captured 
by the o(l) in (08}. 

In a theory with a non-zero mass, the error ei oc due to imperfect particle separation shrinks 
exponentially with the distance between the particles ( §5.5|) . Therefore the total simulated volume 
V should increase logarithmically with 1 / ei oc , and correspondingly the complexity of the algorithm 
scales as poly(log(l/ei oc )). Similarly, by the analysis of £ 14,11 the number of qubits per site scales 
only logarithmically with etmno where 1 — etrunc is an inner product between the exact quantum 
state and the achieved state. Thus, this source of error also contributes only a polylogarithmic 
factor to the overall e scaling. 

By (|139p . the errors resulting from use of a /c th -order Suzuki- Trotter formula with n timesteps 
are esT ~ n~ 2k . est is an operator norm of the difference between the achieved and exact unitary 
transformations. It thus induces a Euclidean distance of O(esT) between quantum states. Hence, 
this error contributes only a factor of 0(e _1 / 2fc ) to the overall e scaling. 

The number of quantum gates used to simulate the strongly coupled theory has scaling in 
1/(A C — Ao) and p that is dominated by adiabatic state preparation ( §4.2.2p . We also estimate 
scaling with n ou t as follows. For two incoming particles with momenta p and — p, the maximum 
number of kinematically allowed outgoing particles is n ou t ~ P- Furthermore, one needs V ~ n out 
to obtain good asymptotic out states separated by a distance of at least ~ l/m-o- For continuum 
behavior, p = rj/a for constant 77 <C 1. Thus, V ~ Hou/- Hence one needs n 2 ^^^ gates to prepare 
the free vacuum and, by (|136p . n 2 ^^ ^ gates to reach the interacting theory adiabatically. (The 
adiabatic turn-off takes no longer than the adiabatic turn-on.) Thus the total scaling in n ou t is 
dominated by preparation of the free vacuum in three-dimensional spacetime, but by adiabatic 
turn-on in two-dimensional spacetime. The results of our resource analysis are summarized in 
Table U 
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K — Ao 


P 


"-out 


d= 1 




p 4+o(l) 




(a c -A J 


d = 2 


/ x n5.U4+o(1) 


p 6+o(l) 


0(n 7 u \ 28 ) 


\\ c -\o ) 



Table 1: The asymptotic scaling of the number of quantum gates needed to simulate scattering in the strong- 
coupling regime in one and two spatial dimensions is polynomial in p, the momentum of the incoming pair 
of particles, A c — Ao, the distance from the critical coupling, and n ou t, the maximum kinematically allowed 
number of outgoing particles. Note that V is kept fixed in calculating the scaling with p. This is justified 
when one is interested in scattering processes with a bounded number of outgoing particles. 

4 Analysis of Algorithm 

Analysis of the algorithm requires quantifying various sources of error. In broad terms, these fall 
into two categories: field-theoretical cutoffs to render the problem finite, and quantum computing 
primitives upon which the algorithm is built. 

Cutoffs are imposed on both space and the field itself. In §4.11 we analyze the effect of dis- 
cretizing and imposing a cutoff on the field, and thereby determine the number of qubits that is 
sufficient to represent the field. The effects of putting space on a lattice are considered later, in 

Sections I4.2H4.4I address quantum computing primitives. In §4.21 we analyze the adiabatic 
preparation of interacting wavepackets. This method induces a phase shift, whose effect on a 
wavepacket is to cause undesirable broadening and propagation. The phase is shown to be pro- 
portional to t/J 2 , where J is the number of adiabatic steps and r is the period of the turn-on. 
The finite period for turn-on causes errors through imperfect adiabaticity. Such so-called 'diabatic' 
errors fall into two classes: particle creation from the vacuum, and the splitting of one particle into 
three. Using the adiabatic theorem, we derive the probability of such events. Two criteria, keeping 
both propagation and diabatic errors small, then determine satisfactory choices of J and r, and 
hence the gate complexity. The physical mass as a function of the coupling features prominently in 
our calculations. For the weak-coupling regime, its form is obtained by perturbation theory ( §5.ip . 
For the strong-coupling regime, we use its known behavior near the phase transition. 

The time evolution during adiabatic turn-on and turn-off, and during the scattering, is imple- 
mented with a Suzuki- Trotter formula. In £ j4.3|, we show that a /c th -order Suzuki- Trotter achieves 
linear scaling in the number of lattice sites, provided that the Hamiltonian is local. This result is 
also used for the phase estimation with which occupation numbers are measured (see §4.4p . 

4.1 Representation by Qubits 

To represent the quantum state of the field, we assign a register of qubits to store 0(x) at each 
lattice point x E Q. Each </>(x) is in principle an unbounded continuous variable. Thus, to represent 
the field at a given site with finitely many qubits, we cut off the field at a maximum magnitude 
(f> max and discretize it in increments of 5$. This requires 

n b = [log 2 (1 + 20 max /^)] (49) 
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qubits per site. In this section we show that one can simulate processes at energy scale E, while 
maintaining 1 — et rU nc inner product with the exact state, with rib logarithmic in I /a, l/e trunc , and 
V. Our analysis is nonperturbative, and thus applies equally to strongly and weakly coupled 4 
theory. 

Let be the state, expressed in the field representation, namely, 

/oo roo 
d(j)v d(f) V V>(0i,---,<MI0i,---,<M, (50) 
-oo J —oo 



and let 



Then 



-'max 



.■'max 



IV'cut) = / <%■■■ / d(j) V 'ip{4>i,...,4> v )\(j) l ,...4> v ) . (51) 



/<Pmax / > <Pmax 
dxf>v / d(p v p((pi,. ■ ■ ,4>v) , (52) 
'0max J— 0max 

where p is the probability distribution 

P {^...,<h) = M<j>i,..., 0v)| 2 • (53) 

In other words, (V'lV'cut) = 1 — Pout, where p out is the probability that at least one of (f>i, . . . , (f)y is 
out of the range [ — (ftm&xi ^max]' By the union bound (Pr(^4 U -B) < Pr(^4) + Pr(.B)), 

(^IV'cut) > 1 - V maxp out (x) , (54) 

where p ut( x ) is the probability that </>(x) is out of the range [— <p max , </> m ax]- 

Let /U0( x ) and cr<f>( x ) denote the mean and standard deviation of </>(x) determined by p. By 
Chebyshev's inequality, choosing ^ max = At$( x ) + C0 V>(x) ensures 

Pout(x)<^. (55) 



Thus, choosing 



max = O |max (^ (x) + /^^(*)) ) (56) 



ensures (V#cut) > 1 - e trU nc- 

Next, we observe the following. 

Proposition 1. Let p and q be Hermitian operators on L 2 (R) obeying the canonical commutation 
relation [q,p] = it. Then the eigenbasis of p is the Fourier transform of the eigenbasis of q. 

Proof. Let f{5) = e'^pe^ 5 . Then /(0) = p and 

±f(8) = e-* s (-i[q,p])e i 1 s = l-, (57) 

therefore 

f(S)=p + Sl, (58) 
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and thus 

frM s = e* 9 * (p + 81) . (59) 
If \p) denotes the p eigenstate with eigenvalue p, then 

pe iqS \p) = e iqS (p + Sl) \p) = {p + 5) e iq6 \p) , (60) 

that is, \p + 5} = e iqS \p) is the eigenstate of p with eigenvalue p + 5. It follows that 

roo roo pea roo 

e iqS dp\p)e~ iqp = dp\p + 5)e- iqp = dp\p)e~ iq{p - &) = e iqS dp\p)e~ iqp . (61) 



Expanding both sides to linear order in 5, we conclude that J^^dp^e tqp is an eigenstate of q 
with eigenvalue q. □ 

By Proposition [TJ the eigenbasis of a rf 7r(x) is the Fourier transform of the eigenbasis of <^(x). 
Thus, discretizing </>(x) in increments of ^(x) is roughly equivalent to the truncation — 7r max < 
7r(x) < 7r max , where 

1 

TTmax = — 77 • (62) 

~ a 0, 



(I' 



0(x) 



By the same argument used to choose (/> max , choosing 



vr max = O (max ^ (x) + a <x) yjJL j j (63) 

ensures fidelity 1 — e trU nc between \ip) and its truncated and discretized version. 

To obtain useful bounds on </> max and 7r max , we must bound /^,( x ), 0"0( x )i M7r(x)j an d ^(x)- To 
this end, we make the following straightforward observation. 

Proposition 2. Let M be a Hermitian operator and let be a quantum state. Then \{\j)\M\\j))\ < 

Proof. For brevity, let (Q) = (i/)\Q\i/)} for any observable Q. The operator (M — (M)l) 2 is positive 
semidefinite. Thus, 

< l{M-{M)lf\ (64) 

= (M 2 -2{M)M + (M) 2 1> (65) 
= (M 2 ) - (M) 2 . (66) 

□ 

Applied to the definitions 

= (^l^(x)lV') , (67) 
^(x) = V^l^x) 2 ^)-^!^)^) 2 , (68) 

^(x) = (V'kCx)!^), (69) 



a 



7r(x) 



v^Kxm-^KxM 2 , (70) 
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Proposition [2] implies that Hqua and c^( x ) are each at most y (tp\(j)(x) 2 \ip) , and //„-( x ) and are 
each at most ('0|7r(x) 2 |'i/'). Thus, by (fSUj) and 



O max 



V 



x£Q V £trunc 



(71) 

TT.nnx O f >» HX , / ' | 77 ! X j 2 1 r ) | . (72) 



V 



x£f2 V ^trunc 



so that, by ([491) and (J62 

n b = 



O ( log ( a max 



V(^Kx)2|V;><V#(y) 2 |V} 



(73) 



etrunc x,yeO 

To establish logarithmic scaling of n b , we need only prove polynomial upper bounds on (ip\(f)(x) 2 \tp) 
and (-0|vr(x) 2 |^). Rather than making a physical estimate of these expectation values, we prove 
simple upper bounds that are probably quite loose. In the adiabatic state preparation described in 
H4.21 the parameters m\ and Ao are varied. The following two propositions cover all the combina- 
tions of parameters used in the adiabatic preparation and subsequent scattering of both strongly 
and weakly coupled wavepackets. 

Proposition 3. Let H be of the form shown in Suppose tuq > and Ao > 0. Let 
state of the field such that (ip\H\vp) < E. Then Vx £ Vt, 

2 W 

<V#(x) 2 |V) < " 



be any 



Itt(x) 



< 



a a rriQ 
2E 



Proof. 



E > 



\H\ 



yen 



trio 



An 



4! 



> 



\**^<k*) 2 \ 



(74) 
(75) 

(76) 
(77) 

(78) 



where the last inequality follows because all of the operators we have dropped are positive semidef- 
inite. This establishes (|74p . Similarly, we can drop all but the 7r(x) term from the right-hand side 
of ([77)1 . leaving 

£><^|a^7r(x)V>, (79) 



which establishes (ITBl. 



□ 



Proposition 4. Let H be of the form shown in |^J). Suppose < and Ao > 0. Let be any 

state of the field such that {ip\H\ip) < E. Then Vx 6 Q, 

(V#(x)V> < — +»/— ^ + 



An 



36m^ 



24 



Xoa d 



E + 



3(V - a d )m* 



2A n 



|vr(x) 5 



3Fm^ 
^A7 



(80) 
(81) 
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Proof. The operator 



I7(x) 



^-0( x ) +4f^( X ) 



(82) 



is sufficiently simple that we can directly calculate its minimal eigenvalue U m i n . If m\ < and 
A > then 

3tti,q 



U u 



2A, 



o 



Thus, for any state 



Hence, recalling (|4]), we obtain 



yen 



> 



-3Vm 4 
2A n 



£ > 



> 



> 



\Hty) 

yen 

yen 
ly^(x) 2 | 



ivr(y) 2 + i(V a 0) 2 (y) + ^0(y) 2 + ^0(y) 4 



^(y) 2 + i(V a ^) 2 (y) 



3Fm£ 



2A n 



(83) 
(84) 

(85) 
(86) 

(87) 



2An 



((871) follows from ((84"]) . ((88]) holds (for any choice of x) because all of the operators we have dropped 
are positive semidefinite. This establishes (|8ip . 

Similarly, dropping positive operators from ((86]) and using ((84]) yield, for any x, 



ml 



(x) 2 + ^>(x) 4 



< [E + 



3(F - q d )mg 
2An" 



(89) 



Applying Proposition [2] with M = </>(x) 2 shows that (■ip\(j)(x) 4: \ip) > (?/>|(/)(x) 2 |V>) 2 . Thus, 



^<V#(x)V> + ^(</Mx) 

Via the quadratic formula, this implies (|8U[) . 



An 



< E + 



3(V - a d )ml 
2Xo 



(90) 
□ 



4.2 Adiabatic Preparation of Interacting Wavepackets 

In this section, we analyze the adiabatic state-preparation procedure. To analyze the error due to 
finite r and J, we consider the process of preparing a single-particle wavepacket. By the analysis 
in £ 15.51 the procedure performs similarly in preparing wavepackets for multiple particles, provided 
the particles are separated by more than the characteristic length 1 /m of the interaction. In what 
follows we use p = |p|, rather than the .D-vector, as will be clear from the context. 

The phase induced by Mj on the momentum-p eigenstate of H(s) (with energy E p (s)) is 

»M = (B, (i±l) + E P (I)) Jl . (91) 
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Taylor expanding E p about s = (j + \)/J yields 



Thus, the total phase induced is 



i(p) 



0{p) 



d 2 E n 



12J3 ds 2 



+ o(,r 



5>(p) 

3=0 



12J 2 







da 



t dE n 



12J 2 ds 



(92) 

(93) 

(94) 
(95) 



where the approximation holds for large J. For a Lorentz-invariant theory, E p (s) must take the 
form 

(96) 



E p (s) = y / p 2 +m 2 ( 



This should be a good approximation for the lattice theory provided the particle momentum satisfies 
p < I /a. Substituting jM]) into ([95]) yields 



0(p) 



dm? 
ds 



24J 2 ^p 2 + m 2 {s) 



(97) 



Next, we consider the effect of this phase shift on a wavepacket centered around momentum 
p. If the wavepacket is narrowly concentrated in momentum, then we can Taylor expand 9{p) to 
second order about p: 

' " ' 08) 



6(p) ~ 6{p) +V-{p-p) + -B-{p-p) 2 , 



where 



V 



B 



89 
dp 

dH_ 

dp 2 



(99) 
(100) 



The phase shift e iT> '^ p ^ induces a translation (in position space) of any wavepacket by a distance 
D. Similarly, the phase shift e l 2 B '(P~P) 2 governs broadening of the wavepacket. From (|99p and (|97p . 
we have 



Tp 



dm 2 
ds 



24 J 2 (p2 + m 2( s ))3/2 



s=l 



s=0 



(101) 



We next determine the complexity by demanding that the propagation length T> be restricted 
to some small constant, and that the probability of diabatic particle creation be small. Together, 
these criteria determine J and r. We can obtain a tighter bound in the perturbative case than in 
the general case, so we treat these separately. 
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symmetric ( 
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Figure 1: The dashed line illustrates schematically the location of a quantum phase transition of 4> 4 theory in 
two and three spacetime dimensions. A and B denote weakly and strongly coupled continuum-like theories, 
respectively. We prepare them adiabatically by following the arrows starting from the massive free theory 
(mj) > 0, Aq = 0). To maintain adiabaticity the path must not cross the quantum phase transition. 



4.2.1 Weak Coupling 



In the perturbative continuum limit a — > 0, m 2 , is negative. For fixed small a, we can adiabatically 
approach a perturbative continuum-like theory by taking the straight-line path depicted in Fig. Q3 
namely, the following parameterization of ()34|) : 

m 2 (s) = (mW) 2 + sA /u, 

A (s) = sX . (102) 

It is particularly efficient to choose fi in accordance with (j!67|) . so that, at first order in Ao, the 
physical mass remains fixed at for all s. In the perturbative regime, this should ensure that 
the path does not cross the quantum phase transition. 

To calculate the variation of physical mass with s, we must go to second order in Ao- The result 

is 



m 2 (s) = (m (1) ) 2 + s 2 m\ 



where ( §5.ip 



[ 0(AS/(mW) a ) , 

0[\\ log (m« fl )) 
I 0{\l/a 2 ) , 



d 
d 
d 



1, 
2, 
3. 



Substituting ([TPS]) into (fTUTI) yields 

Tp 



ITln 



12J 2 



(p 2 + (m( 1 )) 2 +ml) 3/2 



(103) 



(104) 



(105) 



If we are considering a fixed physical process and using successively smaller a to achieve higher 
precision then, by (|104p . it suffices to choose J to scale ad3 



J 



O 

6 
6 



_rp_\ 
XqD) 

arp \ 



d 

d 
d 



1, 

2, 
3. 



(106) 



7 For d — 1 we have simplified the expression assuming p>m 



(i) 
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Note that, for d = 3, J is suppressed by y/a. This is because, as s increases, the (uncancelled) 
two-loop contribution to the physical mass makes the particle very heavy until s is very close to 
one. Hence, the particle propagates slowly, and less backward evolution is required. 

To determine r, we next consider adiabaticity. Let H(s) be any Hamiltonian differentiable 
with respect to s. Let \4>i(s)} be an eigenstate H(s)\4>i(s)) = Ei(s) separated by a non-zero energy 
gap for all s. Let \ipi(t)) be the state obtained by Schrodinger time evolution according to Hit/r) 
with initial condition |^/(0)) = |<^(0)). The diabatic transition amplitude to any other eigenstate 
H(s)\0 k (s)) = E k (s)\0 k (s)) (k ^ I) is [7m 



(Ms)\Mrs)) ~ f da^L^^e^W-wW) (1 + 0(l/r)) 
Jo Ei{u) — hk[a) 

(The integrand is made well-defined by the phase convention (^fc I ^ffi = 0.) Here, 



(107) 



¥>,(«) = / daE^a). (108) 
J o 

In the case that Ei, E k , and (4 > k\^-\ ( l ) i} are s-independent, this integral gives 

(Ms)Mts)) ~ (l - e -(M>) _-^ffi^ )2 (l + 0(l/r)) • (109) 

In the case that these quantities are approximately s-independent, f|109[) should hold as an approx- 
imation. The quadratic dependence on E k — E\ is an adiabatic approximation traditionally used 
in physics. Motivated by applications in quantum computation, mathematicians have developed 
bounds that hold rigorously even with strong s dependence [68]. These more general results have 
a less favorable (cubic) dependence on E k — E\. However, the traditional adiabatic approximation 
appears to be applicable to our analysis. 

In reality, we wish to prepare a wavepacket state, not an eigenstate. However, the wavepacket 
is well separated from other particles and narrowly concentrated in momentum space. Thus, we 
shall approximate it as an eigenstate \4>i(s)). Furthermore, by our choice of path, the energy gap 
is kept constant to first order in the coupling, and thus f)109[) should be a good approximation to 

(unzD. 

Summing the transition amplitudes to some state \<f>k) from the J steps in in our preparation 
process and applying the triangle inequalitjjf] yield the following: 



r 
j=0 



(M/J)\%\M/-0) 



(E k {j/J) - Ei(j/J))< 



(110) 



The j = term in this sum can be evaluated exactly, because it arises from the free theory. At 
j 7^ the theory is no longer exactly solvable. However, one obtains the lowest-order contribution 
to the matrix element (pi, p 2 , P3, P4; s = l|0 4 |vac(l)) in renormalized perturbation theory simply 
by taking the j = expression and replacing ttt-o with the physical mass and Ao with the physical 
coupling. Our adiabatic path (11020 is designed so that the physical mass at s = 1 matches the 



8 The O(J) scaling obtained by the triangle inequality can be confirmed by a more detailed calculation taking into 
account the relative phases of the contributions to the total transition amplitude. 
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bare mass at j = (at least to first order in Ao). Furthermore, the physical coupling differs from 
the bare coupling only by a logarithmically divergent (in a) correction for d = 3 and non-divergent 
corrections for d = 1, 2|£j Thus, we can make the following approximation: 



\{ct> k \Ur))\=o[ J - 



<&(0)|g|*i(0)) 



(i? fe (0)-^(0)) 2 



111) 



Diabatic errors come in two types, creation of particles from the vacuum, and splitting of 
the incoming particles. The matrix element in the numerator of (jllip can correspondingly be 
decomposed as the sum of two contributions. We first consider particle creation from the vacuum, 
approximating \<pj(s)) as |vac(s)). 



By (HMD, 



dH = y a a 

ds ^ 



^ 4 (x) + A /^ 2 (x) 



(112) 



Substituting this into the numerator of (jllip . setting |<^(0)) = |vac(0)), and expanding <p in terms 
of creation and annihilation operators show that the only potentially non-zero transition amplitudes 
are to states |0fc(O)) of two or four particles. The transition amplitude to states of four particles 
arises solely from the cj) A term in The transition amplitude to states of two particles has 

contributions from the </> 4 term and the (j) 2 term in These actually cancel, because of the choice 
of /x defined in ^5.11 (Note that this requires tuning of fx.) At s = 0, the numerator of (jllip is 
therefore 

(pi,P 2 ,P3,P4|^E^ 4 ( X )l VaC (°)) = ^ W +P3+ , 21 ^— ■ (H3) 



4V A /w(pi)w(p 2 )w(p3)w(P4) 



We obtain the probability of excitation due to creation of four particles from the vacuum by 
squaring the amplitude estimated above, and then summing over all allowed combinations of the 
four outgoing momenta. Thus, 



create 



T 2 \ 2 r 
J / H)°Pl+P2+P3+P4,0 



y 

V 2 r 2 (u;(pi) + w(p 2 ) + w(p 3 ) + o;(p 4 )) 4 a;(pi)a;(p2)a;(p3)w(p4) 

Pl:P2,P3,P4Sl 



(114) 



This sum is difficult to evaluate exactly; instead, we shall simply estimate its asymptotic scaling. 
The question is, with which parameter should we consider scaling? There are at least three regimes 
in which classical methods for computing scattering amplitudes break down or are inefficient: strong 
coupling, large numbers of external particles, and high precision. In this section we are considering 
only weak coupling (that is, X/m A ~ D 1), leaving discussion of strong coupling until the next 
section. For an asymptotically large number of external particles, the efficiency of our algorithm 
depends upon strong coupling, for the following reason. A connected Feynman diagram involving 
n external particles must have at least v = 0{n) vertices, so the amplitude for such a process 



By calculations analogous to those in i J5.ll one obtains 



A 



Ao - 
< A - 



3A" 



+ ■ 



for d = 1, 
for d = 2, 

- A + ^4log(m (1) a) + --- , for d = 3. 



+ 
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is suppressed by a factor of ( e £l d ) V , where E is the energy scale of the process. Since E > 
m, many-particle scattering events are exponentially rare at weak coupling, and thus cannot be 
efficiently observed in experiments or simulations. This leaves the high-precision frontier. Recall 
that the perturbation series used in quantum field theory are asymptotic but not convergent. Thus, 
perturbative methods cannot be extended to arbitrarily high precision. 

Hence, in this section we consider the quantum gate complexity of achieving arbitrarily high 
precision. To do so, one chooses a small to obtain small discretization errors, V large to obtain 
better particle separation, r long to improve adiabaticity, and J large enough to limit particle 
propagation as the interaction is turned on. Thus, we wish to know the scaling of -P C reate with a, 
V, t, and J. In this context, we consider m, A, and |pi| to be constants. 

We now estimate the scaling of -P C reate as a — ?■ 0. 



create 



< 



J 2 



V 2 T 2 



E 



pi,P2,P3sr 



V 2 



(w(pi) + W(p 2 ) + w(p 3 ) + w(-pi - p 2 - p 3 )) 4 0;(pi)w(p2)w(p3)w(-pi - P2 - P3) 

3J 2 



V 2 T 2 



E 



pi,P2,P3£r 

|P1|>|P2|,|P3| 



\ 2 
A 



(w(pi) + Cj(p 2 ) + w(p 3 ) + W(-Pl - P2 - p 3 )) 4 Cj(pi)w(p 2 V(p 3 )^(-pi - P2 - P3) 
J 2 



E 



\ 2 



V2t2 „ jH. cr w(pi) 6 w(p 2 )w(p 3 ) 

|Pl|>|P2|,|P3l 



V 2 T 



A 2 

2 pi) ^ 3er ^Piy^R^) 



fj 2 

2 

6 
6 



d d vi 

v 

VJ 2 



VJ 2 
r' 2 a 



d d P 2 



d d P3 



\2 



w(pi) 6 w(p 2 )w(p 3 ) 



d 
d 



1,2, 
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By (TJOSD and gig} , 

^create = O (£) , d= 1,2,3. (116) 

Next, we consider the process in which the time dependence of the 4 term causes a single 
particle to split into three. For this process, the relevant matrix element is 



(P2,P3,P4|yp ^2 a d ^ 4 (x)|pi 



An <5, 



0°P2+P3+P4,P1 



4! 



4V r v / w(pi)a;(p2)w(p 3 )a;(p4) ' 



(117) 



where pi is the momentum of the incoming particle. By our choice of path, the physical mass is 
s-independent to first order in the coupling, and the s dependence of the coupling is only logarith- 
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mically divergent as a — > 0. Thus, by (|110|) . 

p °_ A 0"P2+P3+P4,Pl (11%) 

<PM T 2 V 2 ^ ( W (p 2 ) +a ,(p3) +a ,(p 4 )_ w (p 1 ))4 w (p 1 ) w (p 2 ) w (p 3 ) w ( p4 )- ^ > 

P2,P3,P4£1 

Let us now examine the divergence structure of Psplit as a — >■ 0. In the limit of large volume, 
the sum converges to the following integral: 

— — f d d p2 f d d p 3 — 

t 2 Jr 2 Jr 3 (w(p 2 ) +w(p 3 ) +w(pi - p 2 - P3) - ^(Pi)) 4 w(pi)w(p 2 )w(p3)w(pi - p 2 - P3) ' 

(119) 

If this were divergent as a — > 0, then by approximating the integrand with its value at large |p>2 1 
and |p»3 1 , we would be able to isolate the divergence: 

Psplit ~ I d d P2 I d d P3 1 . (120) 

T MPl) JT Jr (|P2| + |F>3 1 + |P2 +P3|) 4 |P2||P3||P2 +P3| 

However, for d = 1,2,3, this is convergent as a — > 0. Thus, recalling (|106p . we obtain 

P rt -o(£\=l?G)> d=1 ' 2 ' (121) 

We can consider two criteria regarding diabatic particle creation. If our detectors are localized, 
we may be able to tolerate a low constant density of stray particles created during state preparation. 
This background is similar to that encountered in experiments, and may not invalidate conclusions 
from the simulation. Alternatively, one could adopt a strict criterion by demanding that, with 
high probability, not even one stray particle is created in the volume being simulated during state 
preparation. This strict criterion can be quantified by demanding that the adiabatically produced 
state has an inner product of at least 1 — e a d with the exact state. This parameter e a( j is thus 
directly comparable with e tr unc> and the two sources of error can be added. Applying the strict 
criterion, we demand that P sp iit and P C reate each be of order e a d, and obtain 

Strict = 0( — J , d = 1,2,3. (122) 



, e ad, 

Applying the more lenient criterion that P C rea,te/V and P S piit each be of order e a d yields 

Tlenient = O (J-) , d= 1,2,3. (123) 

For a fc th -order Suzuki- Trotter formula, the asymptotic scaling of the total number of gates needed 
for adiabatic state preparation is O ( (Vt) 1+ 2T ) = O ( (Vt /a d ) l+ ^k) ) . Thus, 



V 



2 



2 k 



GSatic = 0\ \ +-) ) • (12 1) 



V \^2k 



^adfabatic ~ O I I - ) j . (12-)) 
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4.2.2 Strong Coupling 



In two and three spacetime dimensions, we can obtain a strongly coupled (that is, nonperturbative) 
field theory by approaching the phase transition ( §2.1.3|) . As in the case of weak coupling, the 
necessary time for adiabatic state preparation depends on various physical parameters of the system 
being simulated, including the momentum of the incoming particles, the volume, the strength of 
the final coupling, the number of spatial dimensions, and the physical mass. To keep the discussion 
concise, we restrict our discussion to the case of ultrarelativistic incoming particles, with coupling 
strength close to the critical value. Under these conditions, the incoming particles can produce 
a shower of many (n out ~ p/m) outgoing particles. Because of the strong coupling, perturbation 
theory is inapplicable, and, even if it could be used, would take exponential computation in the 
number of outgoing particles. 

In the strongly coupled case, we vary the Hamiltonian (|34h with s by keeping the bare mass 
constant at mo and setting the bare coupling to sAo- We choose Ao only slightly below the critical 
value A c , so that at s = 1 the system closely approaches the phase transition, as illustrated in 
Fig. [TJ Examining (|97p suggests that we can estimate phase errors by understanding the behavior 
of m 2 (s) at s = and s = 1, without needing to know exactly what happens in between. From 

cnzD, 



dm 



ds 



& log i -r a 

25.379 A 
16vr 2 



d 
d 



1, 

2, 



(126) 



and, from (fT7|) and (|18 



dm 2 



ds 



-2(A C - A ) , 
-1.26(A C - Ao) ' 26 



d 
d 



1, 
2. 



Thus, (fTPTj) yields 



J 




d = l,2, 



(127) 



(128) 



under the assumption that (A c — Ao) is very small. 

The result (|97h rests on two approximations, a Taylor expansion to second order in (|92h . and 
an approximation of a sum by an integral in (|95[) . The validity conditions for these approximations 
become most stringent at s = 1, where the derivatives of m 2 with respect to s become large. 
Working out the 0(J~ 4 ) term in (|9T|) at s = 1, one finds that it will be much smaller than the 
0(J~ 2 ) term at s = 1 provided 



A,- — Ac 



(129) 



Similarly, higher-order terms in the Taylor expansion are suppressed by additional powers of 
The criterion f)129|) also suffices to justify the approximation of the sum by an inte- 



J(A c -Ao) ' 

gral in (j95|) . 

We must next consider adiabaticity to determine r. In the ultrarelativistic limit, the relevant 
energy gap 7 is ~ This takes its minimum value at s = 1, namely, 



(A c -Ao) 2 

p 1 

(Ac-Aq) 1 



d 



(130) 
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Unlike in the perturbative case, we cannot make a detailed quantitative analysis. However, under 
the condition (|129p . the energy gap 7 changes only slightly in any adiabatic step of the process 
described by ([35]) . ([36% and ([37]) . Thus we apply the traditional adiabatic approximation f)109j) 
and find that each adiabatic step contributes an excitation amplitude of order Depending 
on the relative phases of the excitation amplitudes arising from the J steps, the total excitation 
amplitude could be as high as Indeed, a detailed analysis applying (|109p to (f35|h (f36|h and 



(|37j) suggests that this bound is not overly pessimistic. Thus, to keep the error probability at some 
small constant e aci , we choose 

r~^—. (131) 

We now consider asymptotic scaling with p for fixed Ao- To achieve continuum- like behavior we 
need a <C ±. Thus (fT28|) yields 

J~ T V2 p (d-V/\ d=l,2. (132) 

Substituting (|129p and f| 1 30 [) into (|13ip . we see that we need 

r>p 2 , d = l,2. (133) 

Substituting (|132p and f| 1 30 [) into (|13ip . we see that we also need 

r>p d+1 , d = l,2. (134) 

The scaling r = 0{p d+1 ) for d = 1, 2 suffices to satisfy both conditions (1133P and (I134p . Thus, by 
the results of £ 14.31 the total number of gates scales as 

G strong = O ((yr) 1+ °( V +1+o{1) ) (135) 
= O (V 1+ °(V d+2+ ° (1) ) , (136) 

for d = 1,2. 

Next, we consider asymptotic scaling with (A c — Ao) for fixed p. The J scaling as yfr in (I128P 
automatically satisfies the condition (|129p . Thus, we substitute (|128p into (|13ip . obtaining 

r-{ )^ 5 .04 d = 1 ' (137) 



1 



Hence, using a /c th -order Suzuki- Trotter formula, we obtain 



d = 2. 



Gstrong ~ * 



A c — Ao 
An 



d = 1, 

(138) 

d = 2. 



Note that one could improve this scaling by choosing a more optimized adiabatic state-preparation 
schedule, which slows down as the gap gets smaller. 



32 



4.3 Suzuki- Trotter Formulae for Large Lattices 

It appears that, while scaling with t has been thoroughly studied, little attention has been given 
to scaling of quantum simulation algorithms with the number of lattice sites V. Using a result 
of Suzuki and elementary Lie algebra theory, we derive linear scaling provided the Hamiltonian is 
local. 

For any even k and any pair of Hamiltonians A, B, 

e iAa 1 t/n e iB/3 1 t/n e iAa2t/n e il32Bt/n e iAa r t/n\ n _ e i(A+B)t _j_ Q^2k+1 ^/ n 2fc^ (139) 



where r = l+h k / 2 1 and a\, . . . , a r , /3i, . . . , /? r _i are specially chosen coefficients such that Yfj=i a j = 
1 and J2 r j=i Pj = 1 [SJ. Thus, using the /c th -order Suzuki- Trotter formula (|139p . one can simulate 

/ 2k-\-l \ 

evolution for time t with O ( t~^~ I quantum gates [78]. To determine the V scaling, we use the 
following standard theorem (cf. the Baker-Campbell-Hausdorff formula). 

Theorem 1. Let A and B be elements of a Lie algebra defined over any field of characteristic 0. 
Then e A e B = e c , where C is a formal infinite sum of elements of the Lie algebra generated by A 
and B. 

A and B generate a Lie algebra by commutation and linear combination. Thus, without requir- 
ing any explicit calculation, Theorem [Tj together with (|139p implies 

e iA5 lt /n e iB5 2 t/n _ _ _ e iA8 r t/n\ n = & i{A+B)t + A 2fc+ i* 2fc+1 / n 2k + 0{n^ 2k+1 ^) , (140) 



where A2fe+i is a linear combination of nested commutators. In general, || A2fc+i|| could be as large 
as (max {||A||, H^ll}) 2 ^ 1 . However, by ©, one sees that, for the pair of local Hamiltonians H^, H n , 

/ 2k-\-l 1 \ 

||^2fe+i|| = 0(V), for any fixed k. Thus, one needs only n = O ( t~2^y 21 J . Recalling the O(V) cost 



for simulating each e lH 'i' St or e lH ^ St , one sees that the total number of gates scales as O ( (tV) 1+2k 



Note that this conclusion may be of general interest, as it applies to any lattice Hamiltonian for 
which non-neighboring terms commute. 

In the case of strong coupling, we care not only about how the number of gates scales with V, but 
also about scaling with p. In the presence of high-energy incoming particles, the field can have large 
distortions from its vacuum state. For example, if (^^(x)!^) is large, then local terms in A2/C+1I'!/') 
such as 7r(x)(^(x) 3 |^) can become large. We can obtain a heuristic upper bound on this effect by 
noting that, in the strongly coupled case, rriQ > 0, so each local term in H is a positive operator. 
Thus, if [ij)\H\tp) < E, then the expectation value of each of the local terms is bounded above by 
E. Using E as a simple estimate of the maximum magnitude of a local term, we see that A2fc+i|V , ) ! 
which is a sum of O(V) terms, each of which is of degree 2/c+l in the local terms of H , has magnitude 
at most 0(VE 2k+1 ), or in other words 0(Vp 2kJrl ). Recalling that a scales as a small multiple of 1/p, 
wc sec that Aak+ilVO = 0{Vp 2k+l+d ). Thus, n = 0{p 1+(l + d V 2k t l+l / 2k ). Each timestep requires 
0(V) = 0(Vp d ) gates to implement. Thus, the overall scaling is 0(p d+1+ °^ (tV) 1+ °^) quantum 
gates to simulate the strongly coupled theory at large p. 

4.4 Measurement of Occupation Numbers 

In §2.1.R we saw that the free theory can be decomposed into bosonic modes of definite momentum. 
In this section we analyze the complexity of measuring the occupation of these modes. 
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By (HH), the operator L d a P a p has eigenvalues 0, 1, 2, 3 . . ., which indicate the number of par- 
ticles in momentum mode p. By the method of phase estimation [65], measuring L~ d a P a p reduces 
to simulating e lL da p a P t for various t. By Q, 

L- d a v a p = n p + % + X p , (141) 



where 



*p = y E e*- (x " y) ^^(y). (143) 

x,yeSl 
■ 2d 

*v = ^ E (e ip - (x - y V(x)vr(y)- e ^( x -y)vr(x)0(y)) (144) 

Simulating e li da p a p* with standard Suzuki- Trotter formulae is not very efficient because of the 
0(V 2 ) mutually non-commuting terms in \p- F° r example, the method of [7S] requires 0(V 4 ) 
quantum gates for this task. 

However, there is a simple remedy for this inefficiency. A short calculation shows that 

$p + Up = ( apflp + °-p a - p ) • (145) 

L_ (J , , ,t 



By (|139[) . the problem of simulating time evolution induced by [a p a p + al p a_ p J thus reduces 

to the problem of simulating each of e lIIp * and e** p *. One can efficiently simulate each of e iIlpt and 
e**p* by going to the basis in which it is diagonal, reversibly computing the induced phases, and 
then applying them by phase kickback ( £|2.2.2p . 

The nested commutators appearing in the leading correction to the A; th -order Suzuki- Trotter 
approximation have norm O(V). To see this, first note that by (|142p and (|143p . each power of 



..d 



lip or contributes a coefficient and a summation over V 2 pairs of lattice points. Thus a 
nested commutator at (k + l) th order in n p and ^p consists of a coefficient a -J fc +i ^ times a sum of 
y2(fc+i) nes ^ ec j commutators of </>(x)</>(y) operators and 7r(x)7r(y) operators. Out of these V 2 ( fc+1 ) 
nested commutators, only 0(V k+2 ) are non-zero, because of the product of k delta functions arising 
from ([3]). The sum of 0(V k+2 ) non-zero terms all with coefficient a v k+i ^ yields, by the triangle 
inequality, a total operator norm of 0(V). 

By the above analysis, with a A; th -order Suzuki- Trotter formula, it suffices to use O ^fa^rysE^ 

timesteps. Each step requires 0(V 2 ) gates to simulate, owing to the cost of reversibly computing 
the double sums in n p and <3? p , so the total cost of measuring the combined occupation of modes p 

and — p is O (t 2k y 2 +2fc^. To distinguish eigenvalues separated by a gap 7 using phase estimation, 
we can choose t ~ I/7. All of the eigenvalues of <J> P + n p are separated by gaps of 1/2. Thus, 
t = 0(l). 

This method has the deficiency of failing to distinguish particles with momentum p from those 
with momentum —p. Furthermore, the complexity of this measurement scales quadratically with 
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V. Both of these problems can be corrected by simulating localized detectors. Let Vd C V be the 
region occupied by our simulated detector. As the simplest case, one can consider Vd = olIA ; 
in general, Vd is any spatial translation of this. Then, the corresponding Fourier basis on this 
region is Vd = -P 2 -^ , where Ld = clLd is the length of the detector. A set of operators 

{op,D,Op£)|p £ Td} can be obtained by replacing f2, T, and V with their local counterparts, Ud, 
Td, and Vd, in ©• These operators obey all the commutation relations that one expects of creation 
and annihilation operators. 

Physically, the operator a p D a P ,D can be interpreted as a number operator for the region Vd- 

By the results of §4.51 a p d cl p ,d can be exponentially well approximated by a linear combination 
of cf> and 7r operators with support only within Vd or a distance O(l/mo) of Vd- Thus, by phase 

estimation, we can measure in the eigenbasis of -^TE ( a p D a p, D + a - p D a -p>-°) usrn g 0(V D +2k ) 
quantum gates, independently of V. 

The localized detector ( a P D a p,D + a - p D a -p,£>) detects only particles in the region Vr>. 
Thus it is partially momentum-resolving and partially position-resolving. In accordance with 
the uncertainty principle, the momentum resolution of this detector must be lower than that of 

(a P a p + aL p a_ p ^. This is reflected in the fact that the momentum lattice Td is more coarse- 
grained than r. After a collision, the shower of outgoing particles will be directed outward from the 
collision region. Thus, to resolve the p vs — p ambiguity, it should suffice to surround the collision 
region with a small number of localized detectors. For example, one could simulate 2d detectors 
corresponding to the faces of a (i-dimensional cube surrounding the collision region. 

Number operators with different momenta or within different spatial regions commute ( §4.5p . 
One can thus measure q momentum modes within each of r spatial regions simply by repeating the 

( i+— 

phase estimation procedure for each number operator, with total complexity O I rqV D 2k 
4.5 Localized Wavepackets and Detectors 

In this section we show that the operator a£ defined in (|16|) is quasilocal, namely, it can be expanded 
in the form 

4 = E [/(y " x )^) + 9(y - *My)] , (146) 
y 

where f(y — x) and g(y — x) are each either zero or exponentially small for |y — x| 3> I/tuq. 
This implies that the Hamiltonian defined in f|30[) is quasilocal whenever ip has local support. 
Therefore simulating time evolution according to in order to create a wavepacket as in §3.11 
does not disturb previously created wavepackets provided they are well separated. Furthermore, 
the number of quantum gates needed to simulate iTw, scales only with the number of lattice sites 
in the support of tp, rather than V. 
By (fl6l) and ([9]) one calculates 

4 = I^x) - i £ a d f d (y ~ x)vr(y) , (147) 
yen 



where 



«y-x) = £rt-p*-*>— . (148) 

per yi ^ J 
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In the continuum limit, 

f d d r> p i P'(y~ x ) 

h(y - x) =f / J = . (149) 

For |y — x| ;§> 1/rriQ, this function decays exponentially with characteristic length 1/mo- For 
example, in d = 1 with infinite volume, 

fi{y-x)^ / = -^o(m |y - x|) , (150) 

where -Ko is the modified Bessel function of the second kind, which satisfies 

K (z) ~ i= far z » 1. (151) 

5 Some Field-Theoretical Aspects 

This section describes some quantum field-theoretical details, beginning with the perturbative 
renormalization of the mass in §5.11 

Section 15.21 gives a brief introduction to effective field theory, which is now a well-developed 
formalism underlying our modern understanding of quantum field theory. In its regime of validity, 
typically below a particular energy scale, an effective field theory reproduces the behavior of the 
full theory. Although it consists of infinitely many terms, it can be truncated, with corresponding 
finite and controllable errors. 

In §5.31 we analyse the effect of discretizing the spatial dimensions of the continuum 4 quantum 
field theory. The discretized Lagrangian can be thought of as the leading contribution to an effective 
field theory. From the leading operators left out we can thus infer the scaling of the error associated 
with a non-zero lattice spacing, o. First, in §5.3. H we obtain the general form of the effective field 
theory, including the scaling of the coefficients of different operators. Next, in §5.3.21 we explicitly 
calculate the coefficients of the leading operators in the complete effective field theory at weak 
coupling, matching the full and effective theories at an energy-momentum scale p ~ l/o. The 
effective field theory is demonstrated to consist of three different classes of operators, shown with 
the scaling of their coefficients in Table [2j At strong coupling, the operators and their scaling remain 
the same at the matching scale, although the explicit coefficients are no longer calculable. However, 
the running of the coefficients down to lower energies is determined by their anomalous dimensions, 
which depend on the coupling strength. These anomalous dimensions modify the scaling; at weak 
coupling the modification is small, but at strong coupling it could be larger (though the scaling 
will remain polynomial). 

Section [5.41 addresses finite- volume effects, which should be small, since the interactions in our 
field theory are short-range. This expectation can be confirmed and quantified in the perturbative 
regime. From a technical perspective, the finite volume means that, in the calculations of operator 
coefficients, integrals over loop momenta become (discrete) sums. Consequently the coefficients 
of the operators are modified. The results for the operator are shown in Table [3j Note that 
magnitude of the coefficient is increased for finite length L. 

In §5.51 we assess another effect of a finite volume: state preparation is affected, since asymptotic 
wave-packets are then only approximately free. As before, corrections should be small, since the 
interactions are short-range, and we can confirm this claim in the perturbative regime, this time 
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by calculating the effective potential in the Born approximation. The leading asymptotic behavior 
of the effective potential is given in Table |U 



Class 


Operators 


Scaling of coupling 


I 


4> 2n (n > 3) 


\n a 2n-D 


II 


cpd 2 ^ (I > 2) 


a 21 - 2 


III 


(j > 1, I > 2) 


X j+l a 2j+2l+2-D 



Table 2: Effective field theory operators fall into three classes f i j5.3|) . The general operator in each class is 
shown, with the canonical scaling of its coefficient in D spacetime dimensions. 



D 


Coefficient of 


2 
3 
4 


5 X 3 n 3 
64tt 5 A a 

1287T 5 A a 


^ 3 L 2 J 

lcV2 + 4 |f] 

2(2^ + tt) + |(26V3 + 9tt) Jj] 



Table 3: Wilson coefficient of operator </> 6 /6! in effective field theory for 4 theory in D dimensions, with a 
finite number 2L of lattice sites in each dimension (£ i5.4j) . Corrections to the square-bracketed expressions 
are of order (m 2 a 2 , 1/L 3 ). 



D 


V(r — > oo) 


2 
3 
4 


A 2 1 c ~2mr 
32m 3 ^/irrnr 

A 2 1 c -2mr 
647T 3 ' 2 m (mr) 3 / 2 

A 2 1 -2mr 
1287r 5 / 2 m 3 / 2 r 5 /2 e 



Table 4: Leading asymptotic behavior as r — > oo of effective potential for 4> 4 theory in D dimensions (SJS3J). 
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5.1 Mass Renormalization 



In this section, we calculate the renormalized (or physical) mass of the discretized theory in the 
perturbative regime. First, in £ 15.1.11 we derive the action and the propagator. Next, in §5.1.21 we 
use (a hybrid form of) perturbation theory to obtain the mass to second order in the coupling. 



5.1.1 The Discretized Theory 

In quantum field theory, the Lagrangian formulation is generally more convenient than the equiva- 
lent Hamiltonian formulation. The Lagrangian (density) C((f), d^cft) defines the field theory, and is 
related to the Hamiltonian density H(7r,(fr) by 

H=TT<j>-£, 

where tt = dC/dcft is the conjugate momentum density. For the scalar (ft 4 quantum field theory, 

£ = \d^<ft - - ^Ao^ 4 , (152) 

where fj, = 0, 1, . . . , d in D = d + 1 spacetime dimensions. Then the field and coupling have the 
following mass dimensions: 

M = ^> M=4-U. (153) 
We shall discretize the spatial dimensions, that is, put them on a lattice with spacing a, 

aZ d = {x | Xi/a G Z} . (154) 

The time dimension will be left continuous. Spatial derivatives become (forward and backward) 
difference operators, 

Af/(x) = I + a*) -/(*)) , (155) 
a 

1 



A?/(x) = ~{f{x)- f{x -at)) , 
a 



so that 



V 2 J(x) ee -A<>A{f(x) 
d 



1 a 

i=i 



The free-theory action becomes 

Siree = ~ \ fj dx°dy° £ a 2 <fy(x) (d 2 - V 2 + m 2 ) X)V 0(y) , (157) 



where 

(a 2 - V 2 + m 2 ), (j/ = a- d («9 2 - V 2 + m 2 ) «5 x , y 5(x° - y°) . (158) 
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G(x, y; a) is the inverse of (<9f — V 2 + mo) x , that is, 



J dy° £ a d {d 2 -V 2 a + m 2 ) x y G(y, z; a) = a~^ x , y (-i)5(x° - y°) . 



Using the Fourier transform 

G(x,y;a) 

we obtain the propagator 

G(p;a) 



d_f_ [*/• d d P {x _ y) 



2vr J_ w/a (2tt)< 



G(p; a) , 



W-Eti^sin 2 ^ 



(159) 



(160) 



(161) 



In the limit a — > 0, we recover the familiar propagator of the continuum theory, 

G(p) 



p 2 — TYIq 



(162) 



The Lagrangian of the interacting theory with spatial dimensions put on a lattice is therefore 



(163) 



5.1.2 The Physical Mass 



The analysis of the adiabatic turn-on procedure involves the physical mass. A suitable expression 
for this can be obtained from a hybrid form of perturbation theory, namely partially renormalized 
perturbation theory, in which we use the bare coupling and field but the renormalized mass. From 
(|163p . we have 



£ (0) = ^) j +^-^y-^ 4 -^' 



(164) 



where the mass counterterm is 

S m = Xq/j = tjiq — m 2 . (165) 

The shift in the mass is determined by one-particle irreducible (1PI) diagrams. These are 
diagrams that remain connected after any single line is cut. At first order in Aq, the 1PI insertions 
into the propagator are 



iM(p) 



SI 



+ 



+ 



(166) 



where the second diagram is the counterterm. The calculation of the one-loop diagram, given in 
Appendix B, implies that 



m 2 + Xq/j, , 



8^°g 

. ^ I 

167T 2 a 

(3) 
r 1 
32^ a? 



61 



+ 



+ ... 



for D = 2, 
for D = 3, 
for D = 4, 



(167) 
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where 



r { a ] = 25.379 .. . , 



- . . ,„ - 112.948... . (168) 

Equation (1167P determines how the bare coupling Ao and bare mass tuq are related if the physical 
mass m has a specified value to one-loop order. 

At order Aq, the 1PI amplitude has the additional contributions 



+ 



(169) 



The renormalization condition satisfied at first order in Ao implies that the first two diagrams 
cancel. The calculation of the remaining two-loop diagram (see Appendix B) implies that 



m 



( m (l))2 



A 



+ 



384(m( 1 )) 2 

(mW) 2 + §log(m«a) + 



(m (D)2_^|L4 + 

1536-71-' a* 



for D = 2, 
for D = 3, 
for D = 4, 



(170) 



where mA 1 ) denotes the renormalized (physical) mass at one-loop order, namely, the quantity that 

is kept constant when one follows the path specified by (|167H . (The constant r\ is defined in 
Appendix B.) If, instead, one keeps rriQ constant while turning on Ao, then the expression in terms 
of mo, that is, the result from bare perturbation theory equivalent to the above, is more relevant: 



m 



m2 + :kl 0g 



64 



Q U 



log 



64 



A, 



o 



384m, 



+ 



nr. 



Ao 



" T T6^ a 

(3) 



„(2) 
_1Q__ 

2567T 3 mo 



A | + |log(m a) + --- , 



m '2 I ' Ao I 
m + 32^^ + 



„(3) 



512 



^^log(m a) 



r (3) 

r l ^0 _l_ 
1536vr Y a 7 ~ r 



\2 



for D = 2, 
for I? = 3, 
for D = 4. 



(171) 



5.2 Effective Field Theory 

The formalism of effective field theories (EFTs) is typically used to calculate observables in phys- 
ically relevant theories and hence make predictions. However, the influence of the EFT approach 
extends beyond just providing a tool for tackling otherwise intractable problems: indeed, it has 
profoundly changed our understanding of renormalizability. In this work, the EFT framework is 
applied to determining the scaling of lattice errors. Somewhat similar ideas were employed by 
Symanzik in the construction of improved actions in Euclidean lattice theories [791 f80| [8~T] . 

An effective field theory can be regarded as the low-energy limit of the fundamental theory 
under consideration. An EFT for a full theory is thus somewhat analogous to a Taylor series for a 
function. The canonical example is Fermi theory, in which the four-fermion Hamiltonian is 

n eS = ^(l L Yv L )(u L ^d L ) + h.c. (172) 

The modern interpretation is that this is an effective low-energy theory, in which the W boson has 
been removed as an explicit, dynamical degree of freedom. Pictorially, this corresponds to 
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In other words, terms of order k 2 jMyr have been neglected in the propagator: 



; +<>(-£-). (i7:i) 



tf-M 2 , M 2 , \M W 

An EFT is constructed from only the relevant infrared degrees of freedom and involves an 
expansion in some suitable small parameter. The resulting effective Lagrangian will typically take 
the form 

C cS = E^^ + EEf^ (174) 

n>0 n>l i„ 

where Cj n are dimensionless coefficients, A is the fundamental mass scale below which the EFT is 

In) 

valid, and the local operators (D\J have the same symmetries as the underlying theory. This is 
an infinite series, but the higher the dimension of the operator the more powers of A by which it 
is suppressed. In other words, the lowest-dimensional operators will be the most important ones. 
Thus, in practice, one can truncate the series at some order dictated by the desired accuracy, so 
that one is left with a finite number of operators and hence a finite number of parameters Cj„ to 
determine. If the underlying theory is known and weakly coupled, one may be able to compute the 
parameters. Otherwise, one can take them to be experimental inputs. 



5.3 Effects of Non-zero Lattice Spacing 

In this section, we determine the infinite series of operators comprising the (complete) effective 
field theory whose leading terms are £(°\ The operators not included in fall into three classes: 
operators of the form cf) 2n , Lorentz-violating operators arising solely from discretization effects, and 
Lorentz-violating operators due to discretization and quantum effects. The Wilson coefficients of 
the ignored operators give the error associated with using 

£(0) 

on a spatial lattice to approximate 

the continuum theory. 



5.3.1 The General Effective Theory 

The full (untruncated) effective Lagrangian will have every coupling respecting the (ft — > —<ft sym- 
metry, and so will take the form 

£ eff = £ (0) + + d<p&<\> + + ■ ■ • . (175) 

This can be simplified. First, the chain rule and integration by parts (with boundary terms dropped) 
can be used to write any operator with two derivatives acting on different fields in the form (ft n d 2 <ft. 
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For example, 



tfd^d^ = \d^)d^ -> - l -4?d 2 <p. (176) 

Such an operator can then be simplified via the equation of motion [82, 83J: if this were d 2 (j)+m 2 (j) = 
0, it would imply that the operator (j) 3 d 2 (j) was redundant and could be eliminated entirelvl^l 
Equivalently, one can think of this reduction as making a field redefinition on <j) to eliminate the 
operator: if we write 

C eS = l -(d^f - \m 2 <\? - ^ + Wl 6 + W2 3 a 2 , (177) 
where 77 is a small parameter, then the shift of variables (/>—>•(/> + ?y<72<^ 3 induces 

C e & -»• \{d^f - \m 2 <j> 2 - \^ + vg'i<l> 6 + 0(v 2 ) ■ (178) 

One can then iterate the process, repeatedly shifting variables to remove redundant terms order by 
order in ij. 

Equation (|153j) implies that 

[c] = 6 - 2D , [c] = 8-3D. (179) 

In D = 4 dimensions, [c] = —2 and [c"] = —4. Since the only pertinent dimensionful parameter is 
the lattice spacing, that is, A ~ vr/a, this means that c ~ a 2 and c" ~ a 4 . We see then that, of the 
operators not included in the Lagrangian £(°\ (f> e is more significant than <f> 2n , n > 3, that is, its 
effects are suppressed by fewer powers of pa <C 1. 

In D = 2, 3, the scaling of the coefficients with a is somewhat less obvious, because now the 
coupling A provides another dimensionful parameter (recall that [A] = 4 — D). To obtain the scaling 
of c, one should consider the Feynman diagram that generates the corresponding operator. This 
involves three (p 4 vertices, so 



X 3 a 6 -° . (180) 



(Other diagrams involve higher powers of A and hence are suppressed by higher powers of a.) 
Likewise, the coefficient of (ft 8 will scale as A 4 a 8_D , which means that it is suppressed by a 2 relative 
to the coefficient of <^> 6 . 

In the following subsection, we verify these scalings in the perturbative regime by explicit 
calculation. 

5.3.2 Matching 

We must match the full theory on to the effective theory at a suitable energy scale. What this means 
is that we calculate matrix elements in the full theory and in the effective theory. Comparing these 




10 The situation is modified by certain discretization effects, a subtlety we shall describe later. 
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gives the Wilson coefficients (that is, the coefficients of the terms in the effective Lagrangian). At 
weak coupling, that is, for sufficiently small values of the coupling, this matching can be done in 
ordinary perturbation theory. In our case, the continuum theory and discretized theory correspond 
to the full and effective theories, respectively, and the scale for matching is determined by the 
lattice spacing, a. 



Operators induced purely by Discretization First, consider the matching of the two-point 
function. Taylor expansion of (11561) gives 

d 

-V 2 J(x) = -V 2 /(x)-^-5 l 4 /(x)a 2 + --., (181) 

i=i 

where V 2 denotes the discrete Laplacian and V 2 denotes the continuum Laplacian. Thus, there 
are Lorentz-violating operators induced purely by the discretization. The leading operator of this 
kind is Yli=i 'Pdf^, which we shall denote by a box on a line. Diagrammatically, the matching 
corresponds to 

with the full (effective) theory on the left-hand (right-hand) side of the equation. Letting c/2 be 
the coefficient of the operator Ylt=i ^dfcp = 0$x^ an d expanding the denominator of (|16ip to order 
Px> we then have 



p 2 - m 2 p 2 + a£ p 4 _ m 2 ^ WP * ( p 2 + o2 4 _ m 2)2 



+ *cPx , ; . j\ ^ (183) 



12-^x \f 1 

: + ^ +... (184) 



p2 _ m 2 \2(p 2 — m 2 ) 2 x (p 2 — m 2 ) 2 
a 2 

=> c = --. (185) 

We can write down any operator of this type and calculate its Wilson coefficient (at the matching 
scale) in a similar manner. This matching calculation is independent of the value of the coupling. 
Indeed, it has no direct relation to quantum mechanics: one can think of such operators as arising 
simply because the difference operators in the discretized theory are only approximately equal to 
the derivatives in the continuum theory!^! 
For convenience, we shall use the notation 

X = + ■ + • • • . (186) 



Operators induced by Quantum Effects As mentioned in §5.1.11 quantum effects (which 
correspond to loop diagrams) induce operators of the form cj) 2n , n > 3 in the EFT, but not operators 
of the form </> 2n , these being "protected against" by the <f> — > —<p symmetry. Furthermore, the 
larger n is, the greater the suppression in powers of a. 

"However, their Wilson coefficients will depend on the scale and mix (couple with other coefficients), in accordance 
with renormalization group equations. 
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To obtain the coefficient of the operator </> 6 , we must calculate the Feynman diagram 



A 

X T + permutations , 



(187) 



but it suffices to consider zero external momentum, since this is a non-derivative operator. 
Diagrammatically, matching of the six-point function is equivalent to 



+ 




(188) 



where the diagram on the left-hand side refers to the full theory and those on the right-hand side 
to the effective theory. (Permutations of the loop diagrams are implicitly included.) The result to 
leading order in a (see Appendix C) is that the coefficient c of the term (c/6l)(f> 6 is 



for D = 2, 

= \ -1^ A3 « 3 > for ,0 = 3, 

, -g^(2^/3 + 7r)A 3 a 2 , for£> = 4. 



(189) 



Lorentz- Violating Operators induced by Quantum Effects Recall that we can use the 

equation of motion to eliminate operators of the form (j) n d 2 (j). Note, however, that discretization 
modifies the equation of motion from d 2 4> + m 2 (f) = to d 2 4> + m? cf> — j^a 2 4> 3 d^(f) + ■ ■ ■ = 0. Thus, in 
addition to the Lorentz-violating operators induced purely by the discretization, there are Lorentz- 
violating operators induced by loop effects (combined with the discretization). The leading operator 
of this kind is (j^dfcj) = (j) 3 d*<j). 

The precise coefficients of such operators are technically difficult to calculate even perturbatively 
(since the corresponding external momenta are non-zero). However, one can determine their scalings 
in a, in the same manner as for the other operators: the coefficient of 4> 3 d^(f) scales as X 2 a 8 ~ D . Hence, 
this operator is suppressed relative to the leading operators in the other two classes. 

5.3.3 Strong Coupling 

At strong coupling, perturbation theory is inapplicable and we can no longer calculate the Wilson 
coefficients explicitly. The a dependence of the coefficients at the matching (energy) scale will be 
unchanged. However, the evolution of the coefficients to lower scales affects the scaling, as we now 
describe. 

For simplicity, consider first the perturbative regime in D = 4. When physical processes involve 
disparate energy scales, logarithms of the ratio of those scales generically appear in calculations. 
If the scales are widely separated, such that the logarithm is of the order of the inverse of the 
expansion parameter, the perturbative expansion will no longer be valid. (In our case, the magni- 
tude of log(ma) should be compared with A 2 /(47r).) This problem is dealt with as follows. After 
one matches on to the appropriate EFT at the first energy scale, the (scale-dependent) Wilson 
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coefficients are run down (that is, evolved) to the next energy scale by means of the renormaliza- 
tion group equations (RGEs), the solution of which resums the large logarithms. This process is 
completely analogous to the case of gauge coupling constants, which run and obey renormalization 
group equations. 

The RGEs are characterized by anomalous dimensions, whose effect is to add to (or subtract 
from) the power of the a dependence of the coefficients. In the perturbative regime, the anomalous 
dimensions are small, suppressed by the coupling. At strong coupling, however, the anomalous 
dimensions are incalculable by known methods, and potentially significant. (Nevertheless, the 
known existence of continuum limits in D = 2,3 indicates that the anomalous dimensions of 
suppressed operators will not override their canonical dimensions.) 



5.4 Effects of Finite Volume on EFT 

Consider now a finite length L of each dimension of the spatial lattice, with L lattice sites, so that 
L = La (since periodic boundary conditions are used). Now the d-dimensional integral over loop 
momenta (in weak-coupling calculations) becomes a d-dimensional sum: 



ir/a 
-ir/a 



n/a d d g 



1 2<V 

(2vr) d a d L d 



(190) 



91 Qd 

where each momentum component takes a finite number of values, given by 

2tt L L L 

Qi = — -ni , rii = , h 1, . . . , 0, . . . , 1 . 

aL 2 2 2 



(191) 



The Feynman-diagram calculation proceeds analogously to that in the previous section; the 
only difference is that we now have a Riemann sum, which we know converges to the corresponding 
Riemann integral as L — > oo. The difference between the Riemann sum and corresponding integral 
is given by the Euler-Maclaurin summation formula [84, 85J, which we therefore use to obtain (see 
Appendix D) 



45 
64vr 5 



3„4 



X 6 a 



_ 20 1 „ f 2 2 mV 1 

1 H — + O m 2 a, — — , — 

3 L 2 V L 2 L 3 



for D = 2 . 



(192) 



Using the Euler-Maclaurin formula iteratively for two- and three-dimensional sums, we also obtain 
5 



64vr 5 



A 3 a 3 



10^2 + S^ + OfmV, 
L 2 



m 2 a 2 1 



L 2 L 3 



for D 



(193) 



and 



15 



128vr 5 



2(2^ + tt) + -(26V3 + 9ir)4- + O { m 2 a 2 



1 



m 2 a 2 1 



L 2 L 3 



for D = 4 . 

(194) 



Note that, in all dimensions, a finite L increases the magnitude of the Wilson coefficient. 
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5.5 Effect of Finite Volume on State Preparation 

The procedure for state preparation uses the fact that each particle can be regarded as isolated 
and free when asymptotically separated. For a finite volume, this is an approximation, but the 
corrections should be insignificant, since the interactions are short-range. 

This claim is quantified by the effective potential V(r) created by the interaction. Using the Born 
approximation from nonrelativistic quantum mechanics, one can equate the scattering amplitude 
iM{2 — > 2) in the nonrelativistic limit with (— i times) the Fourier transform of V(r). For identical 
scalar particles, using the symmetric wavefunction of the two-particle system gives two terms, 
namely, 

V(pi-Pf) + V(pi + p f ), (195) 

where, in the center-of- momentum frame, the initial momenta are ±Pi, and the final momenta are 
±P/- 

Consider, then, the QFT amplitude. To lowest order in the coupling A, one simply obtains a 
repulsive contact term, that is, a term proportional to the delta function. Since our interest is 
in the long-distance behavior, we must go to order A 2 . If the incoming (outgoing) momenta for 
2-body — > 2-body scattering are p,p' (k,k'), then the amplitude can be expressed in terms of the 
Mandelstam variables s = (p + p') 2 , t = (k — p) 2 , u = (k 1 — p) 2 . The one-loop amplitude 




(196) 



can be obtained by a straightforward perturbative calculation. In the nonrelativistic limit and 
center-of-momentum frame, the s-channel contribution vanishes, and the t- and u-channels cor- 
respond to the two terms in (I195D . After taking the inverse Fourier transform, we obtain (see 
Appendix E) 



V (2) (r^oo) = 



6 Conclusions 

In this paper, we have established an efficient quantum algorithm for determining scattering am- 
plitudes in a scalar quantum field theory. In particular, the algorithm uses adiabatic turn-on to 
perform the crucial step of preparing interacting wavepacket states. It also discretizes space; a 
detailed analysis of complexity has addressed, among other issues, discretization errors and the 
continuum limit, a fundamentally important aspect of quantum field theory. 

Our quantum algorithm provides exponential speedups over the fastest known classical algo- 
rithms. Specifically, it applies to both weakly and strongly interacting theories, with a run-time 
that is polynomial in the desired precision, as well as the number of particles and their energy (see 
(|48p and Table [1]). In contrast, standard methods in quantum field theory cannot generally be used 
at strong coupling, or beyond a certain precision. 



„ / e~ 2mr + ••• , fori) = 2, 



Q4TT 3 / 2 m (mr) 3 / 2 



e" 2mr + ---, for£> = 3, (197) 



A 2 1 „-2rar 



' 1287T 5 /2 m 3/2 ^572 



e- 2mr + • • • , for D = 4 . 
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We have focussed upon massive scalar <j) theory in spacetime of four and fewer dimensions. 
In future work, we shall extend our results by considering such problems as fermions, gauge sym- 
metries, and massless particles. Our studies pave the way to a quantum algorithm for simulating 
the Standard Model of particle physics. Such an algorithm would demonstrate that, except for 
quantum-gravity effects, the standard quantum circuit model suffices to capture completely the 
computational power of our universe. 



Acknowledgments: We thank Alexey Gorshkov for helpful discussions. This work was supported 
by NSF grant PHY-0803371, DOE grant DE-FG03-92-ER40701, and NSA/ARO grant W911NF- 
09-1-0442. Much of this work was done while S.J. was at the Institute for Quantum Information 
(IQI), Caltech, supported by the Sherman Fairchild Foundation. K.L. was supported in part by 
NSF grant PHY-0854782. He is grateful for the hospitality of the IQI, Caltech, during parts of this 
work. 



47 



Appendices 
A. Notation 



Notation 


Meaning 


n 


Set of all spatial lattice points 


r 


Set of all momentum-space lattice points 


D = d+ 1 


Number of spacetime dimensions 


x,p 


(i-dimensional spatial and momentum vectors 


x,p 


D-dimensional spacetime and energy-momentum vectors 


0(x) 


The field operator at x 


7r(x) 


The operator canonically conjugate to <^>(x) 


H 


Lattice (j) 4 Hamiltonian 


C 


Lagrangian density 


^p? ^p 


Creation and annihilation operators for momentum mode p of the free theory 


m 


Bare mass, not to be confused with m, the physical mass of a particle 


Ao 


Bare coupling, not to be confused with A, the physical coupling 


a 


Lattice spacing 


L 


Length of the lattice 




Total volume of the lattice (V = L d ) 


L 


Number of lattice sites in one spatial dimension (L = L/a) 


V 


Total number of lattice sites (V = V/ a d ) 


r 


Duration of the simulated adiabatic state preparation 




Discrete Laplacian 



Tabic 5: Notation 
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B. Loop Integrals for Mass Renormalization 

At first order in Ao, the 1PI insertions into the propagator are 

-iM(p) = 9 + 
The one-loop diagram gives 



+ 



SI 



,A f°° r /a r /a d D c 



/oo pir/a p 
-oo J—n/a J — 



.Ao a 



2-D r-K 



4 (2vr) d J_ 



I -, t a (2vr) D (g 0)2 _ Yf t=l £ sin 2 {%) - m 2 



d d q 



Zti 4 sm 2 (f) + m 2 a 2 



(198) 

(199) 
(200) 



In D = 2 dimensions, the integral is 

1 



dq 



* ^4 sin 2 (|) + y 2 



yiT^ 2 v 4 + 



(4 + y 2 ) 



y 2 + 



(201) 
(202) 



where i^T(x) is the complete elliptic integral of the first kind. For D = 3 and D = 4, the integrals 
converge at y = 0, so 



/7T />7T 
-7T </ — 7T 



r (c ° + 



where the ellipsis on the right-hand side denotes higher-order terms in y and 

rj 2) = 25.379... , rj 3) = 112.948... . 



(203) 



(204) 



The renormalization condition M(p = m) = implies that i5 m equals the value of the one-loop 
diagram. Thus, 



m 2 + A // , 







, forL> 


= 2, 


n = < 


4 2) 1 , 

16?r 2 a "•" ' 


for D 


= 3, 




r (3) , 
l 32n- s a 1 ~ ' 


for D 


= 4. 



(205) 



At order A 2 ,, the 1PI amplitude has the additional contributions 

■e- ■ 



+ 



+ 



(206) 
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The renormalization condition satisfied at first order in Ao implies that the first two diagrams 
cancel. The remaining two-loop diagram (with external momentum p) gives 



HAo) 5 



ff d D k dP - 
J J (2^(2^ 



3 jojojo 

where a Feynman-parameter integral has been introduced, with 



)° (fc°) 2 - £, £ sin 2 (<f ) - m* (,0)2 _ E . ^ sin 2 (Ml) _ m 2 

(207) 

(208) 



(p0 + fe + g0 )2 _ E . ^ gin 2 ( afr'+fc'+g') ) _ m 2 

iA 2 / ,1 / ,1 / ,1 J , , /■/■ d D g 1 

— / / / dxdydz5(x + y + z-l) J J ^ D ( ^ )D ■ 



D = x 



(*°) 2 -E^ 2 (^) 



D3 ' 



/ 0\2 41 ■ 2 \ 

) -Z^ Sm (-9-) 



+z 



/- ,0 0\2 4 • 2 r a (P % + k l + q l )s 
(P + * + 9 ) - X) 3 sm (— o — ) 



m 



To evaluate the k° and c/° integrals, one can change variables: 



D = ) 8|2 + ^2 + c(p ) 2_ m 2_ a .J-_ sin 2(_)_ J/ J- 



4 . 2 ,ag\ 



where 



+ — (9°+P°) 

x + z 



9° + 



xz; 



xy + xz + yz: 
/3 = x + z , 

xy + xz + yz 



c 



X + z 
xyz 



xy + xz + yz 



Now, 



d/ 2 



2vr 2vr (/3/f + e/ 2 2 -^ 2 ) 3 



*2£ 



2 



(209) 



(210) 

(211) 
(212) 
(213) 
(214) 

(215) 



x 2vr 2vr (pq E + zq E + Ai)z (216) 

= f°<W r * r ^ e -p(0J^+€^ a ) (217) 

= / dpp 2 e- pA (4vrp/3)- 1 / 2 (47rpO _1/2 (218) 
■/ o 

1 

wPF' (219) 
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Thus, we obtain 

-e- 

where 



2 „4-2d 



i\ z a 



24tt (2tt) m 



1 , , , 5{x + y + z-l) 
dx ay az ■ 



o 



y/xy + xz + yz 



d a k I d a q^, (220) 



A = m 2 a 2 



, k' 1 \ . . 2 I't \ , a _• 2 + ^ + 9* 



C V + x ^ 4 sin 2 (— ) + y ^ 4 sin 2 (— ) + z ^ 4 sin 2 



) • (221) 



i=l 



We shall consider this result at the point p = (m, 0), that is, take the renormalization condition to 
be M(p = (m,0)) = 0. 

As a —7- 0, the momentum integral is convergent in D = 4, but becomes singular in D = 2, 3. 
The singular part can be extracted, and the final result is then 



p=(m,0) 



7^ > 



384m 

-^log(ma) + 

ir (3) A 2 
"l A o _i_ 

^ 15367T'' a 2 " ' 



where 



Hence, 



.(3) 







^/xy + xz + yz 



d 6 k 



for L> = 2, 
for D = 3, 
for L> = 4, 

7T 

d 3 o 



(222) 



[A (ma = 0)] 2 



(223) 



m 



[m 



\2 



+ ••• 



384(m( 1 >) 2 

(m^) 2 + §log(mWo) + 



k (m«) 2 



„(3) 



An 



1 _i_ 



1536vr 7 a 



for D = 2, 
for L> = 3, 
for D = 4, 



(224) 



where denotes the renormalized (physical) mass at one-loop order, namely, the quantity that 
is kept constant when one follows the path specified by (11670 . 



C. Loop Integrals for Matching 

To obtain the coefficient of the operator c/> 6 , we must calculate (at zero external momentum 
since this is a non-derivative operator) the Feynman diagram 




+ perms. 



15(-i\y 



oo fir I a 



oo J — w/a 



n/a d D g 



-15U 3 a 6 - D 
1 45 



OO f'TT 



(2vr) d 16 



OO J —77 

— 7T 



/a (27T) 

* d D q 



D * q 2 — m 2 



i\ A a 



K (2tt) d (q 2 - m 2 a 2 ) 3 

d d q - 

(q 2 + m 2 a 2 ) 5 / 2 



(225) 

(226) 
(227) 
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since 



/ 



dp 



2ir ((p ) 2 - .4 2 ) 3 16A 5 
Now, evaluation of the remaining d-dimensional integral gives 



/7T />7T -1 



y2)5/2 



(228) 



(229) 



f 2(27r 3 +37r^ 2 ) 
3j/ 4 (tt 2 +?/ 2 )3/2 



4 (27r 2 j/+(7r 2 +y 2 )^2n 2 +y 2 ArcCot [y^^+^/vr 2 ] ) 

3 ? /3( 7r 2 +y 2)^/27r 2 +2 / 2 
8?r arctan [?r 2 / A /37r 4 +47r 2 j; 2 +j/ 4 



y 2 \J n 2 +y 2 

The power series of these functions around y = are 

^(y) ' 



for d = 1 , 
for d = 2 , 
for d = 3 . 



(230) 



2vr 10v{2 , n t,.2\ 



for d = 1, 
for d = 2, 



(231) 



^ _2(2yjHr) +0(y2)) for d = 3. 



3j/ 2 



Hence, 



\ 7 + perms. 



45 



Z7T 16 

45 A 3 a 3 



(27T) 7 16 



i 45 \3„2 
(2^16^ " 



4 


1 


ma) 4 




27T 


1 


3(mn) 3 




4-7T 


2 


3(ma) 2 



+ 0(m 2 a 2 



3tt 2 



(232) 

for L> = 2, 

for D = 3, (233) 

for £> = 4. 



Diagrammatically, matching of the six-point function is equivalent to 




(234) 



where the diagram on the left-hand side refers to the full theory and those on the right-hand side 
to the effective theory. (Permutations of the loop diagrams are implicitly included.) The coefficient 
c of the term (c/6!)^ 6 is then 



c = 



1 45 

j^y^-X^a^Io^ma) - [a 6 -%-iH] a=0 ), 



(235) 
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since the expression with a = corresponds to the full theory. Thus, to leading order in a, 

-efsAV, forL> = 2, 

-f^A 3 a 3 , for 15 = 3, (236) 

"6i^(2^ + vr)A 3 a 2 , for D = 4. 



Note that the sum of the propagator (|16ip and the series of Lorentz-violating operators that 
begins with <pd^(p is equal to (|162|) . which thus appears in (I225D . Without these operators, (I16ip 
would appear instead, and the resulting six-point function would take the form 



ci 



mo 



+ 



C2 



{ma)' A 



+ c 3 + 



for D = 2 . 



(237) 



In particular, it would have a term of the form a 2 /m 2 , which is divergent in the infrared, (a can 
be regarded as an ultraviolet (UV) regulator and m as an infrared (IR) regulator.) Since the full 
(continuum) theory corresponds to a = 0, this term is absent from the full theory, as (|233h confirms. 
However, the EFT must reproduce the IR behavior of the full theory, and so must have the same 
IR divergences. Therefore, one would not have the correct effective field theory. We see, then, that 
obtaining the correct IR structure depends upon including the Lorentz-violating operators. 

D. Loop Sums for Matching 

The Feynman-diagram calculation proceeds analogously to that in the previous section, with the 
replacement 



id(y,L) 



2 d TT d 2 



L d 



E - E 



i 



2 +n 2 + 



+ n 2 d ) + r 



5/2 



(238) 



where Id(y) = Id(y, oo). This is simply a Riemann sum, which we know converges to the corre- 
sponding Riemann integral as L — > oo. 

The difference between the Riemann sum and corresponding integral is given by the Euler- 
Maclaurin summation formula [84^ [85] . which we can write in the following form. 



Euler-Maclaurin summation formula 

JV-1 / . 

- y f - 



-N 



/(^-^(/M-/(-))+£^)V^ 



(239) 



for m > 1. Here, are the Bernoulli numbers, B m {x) are the Bernoulli polynomials, and {x} 
x — \x\ denotes the fractional part of x. 



Let 



(q 2 + m 2 a 2 ) 5 / 2 



(240) 
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Then, the one-dimensional sum is given by 

/ n 1/27T\2 7T //?7T\3 

S ( T ) /'wL,+°(( T ) 

= Ji (ma 

Hence, 



1 10vr 3 
L2 3(7r 2 + m 2 a 2 ) 7 /2 



2vr\3 



(241) 
(242) 



\ f + perms. 



-i A 3 a 4 

32vr 



_±_ _ _L + 0(mV) _ 2^L (1 + 0(mV)) + o (-L 



(243) 



which implies 



45 
64?r 5 



3„4 



A d a 



, 20 1 ^ / 9 2 "i 2 a 2 1 
1 H — + O m 2 a 2 , 



3 L 2 



1? L 3 



for D = 2. 



(244) 



To calculate the two-dimensional sum, we use (|239p twice. 



I 2 {ma,L) = ~f ^2 



f2m 



— 7T — 7T — 7T 

T < \ 207r3 ? 1 ^ f M 

3L 2 7 (<? 2 + vr 2 + m 2 a 2 ) 7 / 2 V.L 3 / 

— 7T 

Evaluating the second term (exactly) and expanding the result around ma = 0, we obtain 



(245) 

(246) 
(247) 



=p f + perms. 



• 45 x3 3 

64tt 2 



2vr W2 , n , 2 2 , 4 3^2 1 „ , ,„ , >, 



3(mo) 3 97T 3 
1 



9vr 3 L2 



+0 



1? 



(l + 0(mV)) 



(248) 



which implies 



64vr 5 



A 3 a 3 



z- 43\/2 ^ / 2 2 m 2 o 2 1 A 
10\/2 + + O m 2 a 2 , , — 

L 2 V L 2 L 3 / 



for L> = 3 . (249) 
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Finally, 
h{ma,L) 



2tt 
L 



—-1 

2 x 



/ ( ~r, Q2, 93 ) dq 2 dq 3 + - ( -Z- J / / ?3 I g 2 , tt ] dg 2 + O 



2ni 
L 



3 vl; 



2vri 
L 



77 TV 7V 



f{qi,Q2,q3)d- i q + - [ -z- 



2 V L 



f<i,,{«i-<l-2)d<li<lq-2 ( 



-7T — 7T — 7T 



h(ma) 



lOvn 



7T 7T 

// 



/ 2 I 2 i 2 i 2 o\7/2 dqidq2 + ( "R 

Wl + 3 2 + ^ + 771 a ) \L 3 



1 

Z 3 " 

(250) 
(251) 

(252) 



Evaluating the second term (exactly) and expanding the result around ma = 0, we obtain 




+ perms. 



45 



3„2 



1287T 3 



4tt 2(2^/3 + tt) 



3(ma) 2 



3vr 2 



+ O(mV) 



27vr 2 L 2 VL 3 



,(253) 



which implies 
15 



3„2 



1287T 5 



A^a 



2(2^ + tt) + ^(26^/3 + 9^-1 + (m 2 a 2 , ^ 



for D = 4 . 

(254) 



Note that, in all dimensions, a finite L increases the magnitude of the Wilson coefficient. 
E. Integrals for Effective Potential 

If the incoming (outgoing) momenta for 2-body — > 2-body scattering are p,p' (k,k'), then the 
amplitude can be expressed in terms of the Mandelstam variables s = (p + p') 2 , t = (k — p) 2 , 
u = (A/ — p) 2 . A straightforward textbook calculation (in renormalized perturbation theory) of the 
one-loop amplitude 



iM 




+ 




+ 




+ 





+ 



(255) 



gives 



iM = 



. . ,A 2 r(2-f) r\ 



2 (4 



1 



1 



(m 2 — x(l — x)s) 2 D / 2 (m 2 — x(l — x)4m 2 ) 2 D l 2 
1 1 



+ 



( m 2 _ x (l _ x ) t ) 2 ~D/2 ( m 2)2-£>/2 
1 1 



( m 2 _ x (l _ x ) u )2-D/2 ( m 2)2-D/2 



(256) 
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(Here, the following renormalization condition has been used: iAi = —iX at s = Am 2 , t = u = 0. 
This corresponds to defining A as the magnitude of the amplitude at zero ci-momentum.) For 
example, in D = 4 dimensions, 



iM 



-iX — i 



X' 



327T 2 Jo 



dx 



log 



m 2 — x(l — x)s 
m 2 — x(l — x)Am? 



+ log 



m 



x(l — x)t 



+ log 



m 2 — x(l — x)u s 



In the nonrelativistic limit, p = (m, p), etc., so that s = 4m — |p + p'\ , t = — |k — p| , 
u = — |k' — p| 2 . Thus, in the center-of-momentum frame, the s-channel contribution (namely, the 
first line of the integrand of (|256p ) vanishes, and the t- and u-channels correspond to the two terms 
in (|195p . One must also divide by (y / 2m) 4 to account for the difference in the normalization of 
states. 

Hence, the potential is V(x) = V^\x.) + V^(x), where 



A 



4m 2 



(258) 



and 

y( 2 )(x) 



\2 TO _ D\ f l f -oo AD-\ n 



4m 2 (4tt)^/ 2 



dx 



o 



(27T) 



D-l 



(m 2 + x(l - x)q 2 ) 2 -^/ 2 (m 2 ) 2 - D / 2 



In D = 2 dimensions, 

V {2 \r > 0) 



A 2 



dx 



dq 



iqr 



16nm 2 7 7-oo 2vr (m 2 + x(l — x)q 2 ) 

y2 pi e -mr/y/x(l-x) 

dx- 



32vrm 3 J ^(1 - x) 



A 2 



1 



32m 3 yfwmr 



^—2mr _|_ 



as r — > oo . 



(259) 

(260) 
(261) 
(262) 



In the last line, the asymptotic evaluation of the integral was obtained by Laplace's method. 
In D = 3 dimensions, 



F (2) (r > 0) 



A 2 



1 



3W (2vr) 2 J dX L dqq L d9 



a iqr cos f 



y/m 2 + x(l — x)q 2 



A 2 



dx- 



-mr/^J x(l— x) 



6An 2 m 2 r J ^J x {\ _ x ) 

^ 2 \ c -2mr _|_ . . . 

64-7r 3 / 2 m (mr) 3 / 2 



as r — > oo . 



In D = 4 dimensions. 

A 2 



V (2) (r > 0) 



1 



OO fl 

2 



64A 2 (2.) 3 J dq * L dU eiqm J *<> h dx bg < 1 + 



2tt 



x(l — x) 



m. 



A 2 



647r 2 m 2 (2Tr) 2 ir 



777T 



(263) 

(264) 
(265) 

(266) 
(267) 



(257) 
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To evaluate the integral with respect to q, complete the contour in the upper half-plane. There is 
a branch cut on the imaginary axis from q = 2im to q = ioo (since the logarithm has a branch cut 
where its argument is negative). 
Let 



f(p) ee / dxlog[l-^-^p 2 ) . (268) 



i 

o ~ V rn' 



Then 



oo 

oo 



poo 

2i dqqe-^lm[f{q 2 - ie)] , (270) 

J2m 



where q = —iq. 

For a fixed q 2 , there are contributions to Im/ when x is in the range x_ < x < x+, where 



1 1 / 4m 2 

= 2 ± aV 1 ~"F _ ■ <271) 



Hence, 



Im[f(q z - ie)} = it j ' dx = vr Jl - ^ , (272) 



and 



■ <!< ""' / '/■<' log j I + ^/V l = -2tt» /°° dqqe-* r Jl - . (273) 



2m V 9 



Substituting this result into (|26T[> . we obtain 



"- ,2) <">»> = -jd^C^'T^f (274) 

A " 1 -,e- 2mr + • • • , as r -»• oo . (275) 



128vr 5 /2 m 3/2 r 5/2 
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